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Abstract 


We prove the existence of diffusing solutions in the motion of a charged particle in the presence of an ABC 
magnetic field. The equations of motion are modeled by a 3DOF Hamiltonian system depending on two para¬ 
meters. For small values of these parameters, we obtain a normally hyperbolic invariant manifold and we apply 
the so-called geometric methods for a priori unstable systems developed by A. Delshams, R. de la Llave, and 
T.M. Seara. We characterize explicitly sufficient conditions for the existence of a transition chain of invariant tori 
having heteroclinic connections, thus obtaining global instability (Arnold diffusion). We also check the obtained 
conditions in a computer assisted proof. ABC magnetic fields are the simplest force-free type solutions of the 
magnetohydrodynamics equations with periodic boundary conditions, so our results are of potential interest in the 
study of the motion of plasma charged particles in a tokamak. 
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1 Introduction 


The study of the motion of a charged particle in a magnetic field is a classical subject in several areas of physics, such 
as condensed matter theory, accelerator physics, magnetobiology and plasma physics. The equation of motion of a 
(non-relativistic) unit-mass, unit-charge particle at the position g S in the presence of a magnetic field B is given 
by the Newton-Lorentz law 

q = qx B(q), (1) 

where the dot over q denotes, as usual, the time derivative, and x stands for the standard vector product in 

An important observation is that Eq. ([Hi can be written equivalently in a Hamiltonian way whenever there is a 
globally defined vector potential A such that B = curl A. If this is the case, the Hamiltonian function is 

H{q,p) = ^{p-A{q)f. 

In this paper we are interested in the motion of charges in ABC magnetic fields. These fields arise in the theory 
of magnetic dynamos (see ll29l and references therein) and were introduced independently by Arnold Q and Chil¬ 
dress lITOl in the 1960’s. The well-known family of ABC magnetic fields depends on three real parameters. A, B and 
C, and reads in Cartesian coordinates q = (x, y, z) as 

Babc = {A sin z + C cos y, B sin x + A cos z, C sin y + B cos x). (2) 


ABC magnetic fields are stationary solutions of the magnetohydrodynamics equations of force-free type, thus 
implying that the field exerts no force on the current distribution generating it. Indeed, it is straightforward to check 
that Babc is divergence-free and force-free because curl Babc = ^ABC^ and so the ABC field admits the globally 
defined vector potential Aabc = ^ABC- ABC magnetic fields are minimizers of fhe energy functional f B^ acting 
on the space of divergence-free fields of fixed helicity. 

Since the dependence of the ABC magnetic field and its vector potential with the variables x, y, z is 27r-periodic, 
it is customary to consider that these fields are defined in the 3-torus = R^/(27rZ)^ so that {x,y,z) € T^. By 
rescaling and reordering the space variables and the time, all the non-trivial cases can be reduced to A = 1 > B > 
C > 0, so we shall assume it in what follows. The Newton-Lorentz equation of motion (dll for the ABC magnetic 
field can be described as a 3DOF Hamiltonian system defined in the phase space x 9 {x, y, z,px,Py,Pz) by 
the Hamiltonian function: 
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H = - {px — C cos y — sin zf -\— {py — B sin x — cos zf — {pz — B cos x — C sin y) . 


(3) 


Force free fields are very important in applications and model diverse physical systems, as stellar atmospheres ||71, 
the solar corona |[24l and relaxed states of toroidal plasmas ll^ . Moreover, the motion of a charge in an ABC field can 
be interpreted as a model for the motion of plasma charged particles in a tokamak. A wide examination of system (IH 
was recently presented in [361, proving the existence of confinement regions of charges near some magnetic lines and 
also that the problem gives rise to non-integrability and chaotic motions. In this study we go one step further and 
obtain global instability, i.e. Arnold diffusion. 

Characterizing global instabilities in Hamiltonian systems is a relevant problem that has called the attention of 
mathematicians, physicists and engineers. For example, in the context of beam physics, designers of accelerators 
or plasma confinement devices are interested in the characterization of these instabilities in order to avoid them as 
much as possible (e.g. in the confinement of hot plasmas for fusion power generation, diffusion is a very relevant 
phenomenon because of the harmful plasma-wall interaction). Global instability deals with the question of whether 
Hamiltonian perturbations of a regular integrable system accumulate over time, giving rise to a long term effect, or 
whether they average out. This problem was first formulated by Arnold. Indeed, in the celebrated paper |[T1, Arnold 
constructed a concrete example, suitably and cleverly chosen, such that some trajectories can jump around KAM tori 
thus obtaining diffusion (after |[ll this problem is known as Arnold diffusion). These diffusing orbits were constructed 
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using a mechanism of transition chains. It consists in obtaining heteroclinic intersections between the stable and 
unstable manifolds of a sequence of whiskered invariant tori. 

In the last decades there has been a significant advance in the understanding of diffusion and, following 0, 
the studies are classified in fwo differenf groups: fhe a priori unstable case and fhe a priori stable case. Arnold 
diffusion in a priori unsfable systems (where fhe unperfurbed system has hyperbolic properfies of some kind) has 
been approached using geomefric mefhods in ifT^ [T^ [T^ l25l . fhe separafrix map in l43l . topological mefhods 
in and variafional mefhods in ||3l[IIl. A combination of fopological and geomefric mefhods has been recenfly 

presenfed in ll28l . The more difficulf case of a priori sfable systems (where fhe unperfurbed system is foliafed by 
Lagrangian invarianf fori) is less undersfood, buf significanf advances have been presented along fhe lasf few years 

in 01 [321 [32 [Ml [3l [Ml 1171. 

Our aim in fhis paper is to prove fhe existence of Arnold diffusion in fhe dynamics of a charged parficle in an 
ABC magnetic field, which is modeled by fhe Hamiltonian sysfem ([ 2 - If H = C = 0, we obfain an infegrable 
Hamiltonian system Hq having a normally hyperbolic invariant manifold (NHIM) Aq foliated by whiskered invariant 
tori (see details in Section [2Jl) . Then, the problem considered in this paper falls into the a priori unstable setting. It 
is worth mentioning that one of the main difficulties in the study of a priori unstable systems was the so-called large 
gap problem (see llTOl f. This problem arises because a generic perturbation of size e creates gaps at most of size ^/£ 
between the persisting primary KAM tori and, in principle, only orbits separated an amount e could be connected by 
heteroclinic connections between invariant tori. This issue was solved in the previously mentioned references, using 
different tools for the study of Arnold diffusion. We observe that recent mechanisms of diffusion have been proposed 
in order to avoid big gaps using very little information of the dynamics restricted to the NHIM (see Ell)- Here, 
we follow the geometric methods developed in ifT^ [Ml in order to prove the existence of Arnold diffusion in the 
Hamiltonian ([2 for small values of B and C. Concretely, we prove the following theorem, which establishes suffi¬ 
cient conditions for the existence of a transition chain between whiskered invariant tori, thus producing large unstable 
motions in the perturbed system: 

Main Theorem (informal statement). Let us consider the Hamiltonian ([2 with B = eB / 0 and C = eC / 0, and 
a non-empty setX = [oi, 6 i] x [ 02 , b 2 ] for given (positive) values ofai, bi. Then, under some explicit non-degeneracy 
and transversality conditions, if\£\ is small enough, the ABC system exhibits Arnold diffusion in X, i.e. there exists a 
trajectory of ([2 connecting two arbitrary values of {px,Py) in the interior ofX. 

A precise statement of this theorem is given in Theorem 12.21 (Section [ 2 , after a detailed discussion of the un¬ 
perturbed ABC system. Moreover, we implement the non-degeneracy and transversality conditions included in the 
Main Theorem in a computer assisted proof (CAP) in Section [2 As a consequence, we obtain an open set of initial 
conditions in phase space where we can construct a transition chain. For example, we obtain the following result 
which serves as an illustration: 

Corollary 1.1. Let us consider Hamiltonian ([2 with H = 10 and C = 0.1. Then, the non-degeneracy and transvers¬ 
ality conditions of the Main Theorem hold in the set X = [0.1,0.9] x [0.5,0.9]. Therefore, for ]e] small enough, there 
exists a trajectory of ([2 connecting two arbitrary values of (px,Py) in (0.1,0.9) x (0.5,0.9). 

We remark that the choice B = 10 and C = 0.1 is arbitrary. Analogous results can be obtained for any other 
choice of parameters. The computational cost to verify the hypotheses for a fixed sef X increases when fhe difference 
befween B and C is reduced. If is worth menfioning fhaf if we fake “narrow” sefs of fhe form X = [oi, oi -|- 5] x [ 02 , 62 ] 
or X = [oi, 61 ] X [ 02,02 -h (5], wifh <5 small, fhen fhe compufafional cosf of fhe CAP is reduced significanfly. In fhis 
case, if is also possible fo check fhe conditions for open sefs of paramefers B and C. We have produced analogous 
resulfs fo Corollary 1 1.1 l and we have nof found obsfrucfions fo diffusion in any case. 

To fhe besf of our knowledge, fhe Main Theorem and Corollary 1 1.1 1 are fhe firsl rigorous resulfs on fhe existence 
of diffusing orbifs in fhe motion of charges in magnetic fields, even fhough physicisfs have been aware of fhis phe¬ 
nomenon for a long time (cf. 1451 [^ i and fhe effecl is somefimes known as driff mofion in fhe physics liferafure. Of 
course, we wanf fo menfion ofher significanf problems where Arnold diffusion have been characferized. In particular. 
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we can find remarkable contributions in the context of celestial mechanics: diffusion along mean motion resonances 
in the restricted planar three-body problem |[2^ : instability mechanism in a special configuration of the 5-body prob¬ 
lem 13911481; transition chains of invariant tori around the point L 2 in the elliptic three body problem as a perturbation 
of the circular problem IT^ . improved recently in 0; instability around the point Li in the circular spatial restric¬ 
ted three-body problem, focusing on homoclinic trajectories lT3l : instability in the elliptic restricted problem close 
to the parabolic orbits of the Kepler problem between the comet and the Sun 1^ . We observe that some parts of 
the arguments in l6l [121 [El |22l are non-rigorous, but are strongly backed by convincing numerical computations. 
It is also worth mentioning the example discussed in ll20l . where the geometric mechanism for diffusion introduced 
in l[T^ is illustrated in a representative model. The model simplifies some of the hypotheses, thus saving a significant 
amount of computations, so they can present the geometric mechanism of diffusion in a clear understandable way. In 
the system Q studied in this paper, some of these simplifications cannot be used and we must perform some ad hoc 
analysis and specific computations. The reader interested in numerical studies is referred to |[30l. 

The mechanisms governing Arnold diffusion are very complex and there are still many questions to answer and 
many aspects to understand. As is posed in ll27l . it is relevant to detect, combine, and compare different mechanisms 
of diffusion displayed by concrete systems. In this way, Hamiltonian Q can be an ideal framework to apply and 
compare different approaches and methods in the literature (e.g. topological methods, variational techniques, use of 
multiple scattering maps, etc). On the one hand, the ABC system is complicated enough to contain all the difficulties 
that are present in a general a priori unstable problem. On the other hand, the ABC system is explicit and simple 
enough to perform analytic computations. Moreover, it is a problem that appears in a natural way in physics. 

The proof of the Main Theorem consists in combining the internal dynamics on the NHIM with its outer (asymp¬ 
totic) dynamics, which is modeled by the scattering map |[T7l . The procedure is divided in the following steps: 

Characterization of the NHIM: The first step is to characterize the perturbed NHIM Ag and its stable and unstable 
manifolds (we summarize some basic concepts in Section ISTI) . We pay special attention to describe explicitly 
the geometric procedure that allows us to parameterize the NHIM in a natural way, thus obtaining a suitable 
symplectic structure on the NHIM (see Section 13.21) . The construction presented has special interest since 
we give explicit formulas to use the deformation theory introduced in |[T7l . To this end, we have to compute 
perturbatively a symplectic frame associated to the manifold. Explicit computations for the ABC system are 
detailed in Section [33] 

Invariant tori on the NHIM: To study the inner dynamics on the NHIM, where the so-called big gaps are present, 
we perform averaging theory of the vector field restricted to the manifold. After choosing a suitable paramet¬ 
erization in the previous step, we follow |[T6l [El mutatis mutandis in Section 13.4.11 Explicit computations for 
the ABC system are detailed in Section 13.4.21 In Proposition 13.9! we obtain an approximation of the level sets 
that characterize the invariant objects inside the NHIM. In particular, we find a set of whiskered invariant tori 
(primary and secondary) covering Ag except for a set of measure 0(e^3)_ 

Scattering map: In Section |4| we describe the outer dynamics associated to our problem. Eor the sake of complete¬ 
ness, in Section l4T] we summarize the construction of the Melnikov potential that characterizes the intersections 
of the stable and unstable manifolds associated to the NHIM (cf. lidTl l. In Section l43] we compute the scattering 
map for the ABC system. 

Combination of inner and outer dynamics: The combination of both dynamics, obtaining explicit transversality 
conditions for the existence of diffusion, is performed in Section jS] We remark that, since the unperturbed 
scattering map has a so-called phase shift, there is an additional term in the transversality conditions that is 
not present in |[T6l [TSl . In the domain where the conditions are satisfied, we construct a sequence of 

whiskered tori satisfying VEji iti that is, we construct a transition chain along Ag. 

We remark again that the hypotheses in the Main Theorem are explicit and involve a series of standard, but 
cumbersome, computations. Eirst, we evaluate some integrals that depend on {px,Py) as parameters. We solve a 
one-dimensional non-linear equation that depends on these integrals. We approximate the derivatives with respect to 
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parameters of the previous solution. Finally, we evaluate several complicated formulas that depend on the previous 
objects. In Section^we rigorously perform these computations with the help of a computer. 


2 Setting of the problem and statement of the main theorem 

In this Section we present a detailed description of the geometry of our problem and state a precise version of the 
Main Theorem. More precisely, in Section 12.11 we fully describe the motion of the unperturbed Hamiltonian system 
(Eq. (I3ll with B = C = 0), thus characterizing a normally hyperbolic invariant manifold with coincident stable and 
unstable invariant manifolds. Then, in Section [2^ we provide explicit sufficient conditions for the existence of Arnold 
diffusion in the perturbed problem (Eq. ([3]) with B = eB, C = eC). 


2.1 Geometric features of the unperturbed problem 

For B = C = ^, the ABC magnetic field has fhe simple expression 

Babc = (sin 2 ;,cosz, 0 ), 


which implies fhaf fhe field is linear on each toroidal surface z = zq, periodic or quasi-periodic depending on fhe 
value of tan zq. Concerning fhe equafions of mofion, fhe Hamiltonian funclion in Eq. (I3]l is given by 


Ho = ^iPx- sin zf + ^{Py- cos zf + , 


The system of ODEs associated to (IHl is 


( 4 ) 


X =Px-sinz, Px = 0, 

y = Py - cos z , Py = 0, 

z =pz, Pz = PxCosz - Py sin z , 

so px and py are consfants of fhe mofion. There is no loss of generalify in faking positive values of px and py, so 
we shall assume if fhroughouf fhe paper. In addifion, we observe fhaf fhe sysfem {z,pz) is pendulum-like and has an 
effeclive potential 

V(z) := —px sin z — Py cos z . 

Notice fhaf fhis sysfem has a hyperbolic equilibrium af fhe poinf 

z* := arctan — + vr, p* = 0 , 

Py 

and, since Px> ^ and py > 0, we have fhe identifies 


sin 2* 


-Px 



cos z* 



We denofe fhe posifive eigenvalue of fhe linearized equation af fhe hyperbolic equilibrium as 

>'■= {pl+plf^^, ( 5 ) 

which allows us fo write fhe consfanfs of fhe mofion as px = sin a, and py = cos a, wifh a = aicinn{px/ Py) S 
(0, vr/2). Wifh fhis nofafion, fhe pendulum-like equation in fhe variables {z,pz) reads as 

z = Px COS z — Py sin z = sin(Q; — z ), 
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thus obtaining that there is a homoclinic orbit connecting the equilibrium point given by 

2A 

z^{t) = 4arctane'^* + z*, P^t) = — . (6) 

cosh(At) 

It is straightforward to check that z^{t) —)• z* and (t) p* = 0, exponentially with exponent A, as i ±oo. There 

is a second homoclinic trajectory connecting the equilibrium point given by z^{t) = —z^{t) + 2a and p^z{t) = —p^{t), 

but it will not be used in what follows. 

The previous computations show that the Hamiltonian system Hq has a 4-dimensional normally hyperbolic in¬ 
variant manifold 

Ao := {{q,p) G X : z = z*,p^= p*} , 

which is foliated by 2-dimensional invariant tori Tp^^py obtained by fixing p^ and py, i.e. Aq = IJp^ py %>x,Py- A direct 
computation shows that the dynamics on each invariant torus Tp^^py is linear with frequency vector w = {i 0 i,i 02 ) 
given by 

:= Px - sin(z*) = p^(l -h {pi , (7) 

UJ2 ■■= Py - C0S{Z*) = Py{l + {pI+pI)~^^^) . ( 8 ) 

The stable and unstable manifolds of Aq are 5-dimensional invariant sets defined by 

H"^(Ao) = H""(Ao) = {(g,p) G x : z = z°(t),p, = p°(r),T G M} , 


so the set Vk®(Ao) (or Vk'^(Ao)) is the union of the stable (unstable) manifolds of the invariant tori Tp^,py, i.e. 

H^^(Ao)= U W^{Tp^,py) = W'^{K^)= U W'^{Tp^,py). 

Px iPy Px iPy 


In order to work with the invariant torus Tp^^py and its whiskers W^{Tp^^py) = W'^{Tp^^py), we introduce appro¬ 
priate parameterizations. Indeed, Tp^^py C Aq can be parameterized as 


u* = u*{x,y) = {x,y,z*,p^,py,pl), 

where p^ and py are fixed and {x,y) G T^. Moreover, the stable manifold W^{Tp^^py) is given by the set of points of 
the form 

= «°(r,x,y) = {x + Fi{T),y + F 2 {T),z^{T),py,,py,pl{T)) , (9) 

where r G M, (x, y) G T^, the functions z^ and p° are defined in Eq. dh]), and 

Fi(r) := sin(z*)r — / sin(z°((T))d(T , F 2 (r) := cos(z*)r — / cos( 2 °(iT))(iiT . 

Jo Jo 

Finally, we introduce some notation that will be useful in Section |4~T] If (j)^ is the flow of the Hamiltonian system F[q 
and we consider points u* G Aq and G H7®(Ao) = IE“(Ao), then 

= {X +UJlt,y + UJ2t,Z*,Py;,Py,pl), ( 10 ) 

= {x + Fi{t + t) +UJlt,y + F 2 {t + t) +UJ 2 t,Z°{T + t),Pa;,Py,pl{T + t)), (II) 

We observe that the functions Fi and F 2 depend on the constants px and py through z* and z^^, but we omit 
this dependence in order to avoid cumbersome notation. After straightforward computations we obtain the following 
explicit formulas 


Fi{r) 

F2{t) 


^2(tanh(AT) — 1) ^ 2^ . ^ ^2sech(Ar) 

/ 2(tanh(AT) — 1) 2 A ^ /2sech(Ar) 

(-A-+(-^ 



COS z* , 
sin 2 ;* , 


( 12 ) 
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where the constant A is defined in Eq. ([51). These functions allow us to compute the phase shift of any trajectory when 
traveling along W^{Tp^^py). Indeed, the phase-shift is defined by the limits 

x+ := limt^oo Fi{t + t), X- := limt^_oo -Fi {r + t), 

y+ := limt^oo F2{t + t), y_ := limt^_oo F2{t + t ), 


which can be explicitly computed and do not depend on r, that is 

^ O (FPx - Py) ^ o (PxFPy) 

{pI + ’ {pI + pI)^^‘^ ' 


(13) 


Observe that the limits and x_ are different, which means that any point in the homoclinic orbit approaches 
different points of the same invariant torus if we consider the limit in the future and in the past. This is the reason 
why the terminology phase-shift is used for this phenomenon, see e.g. |[5l[T5l|22l. As we will show in Section [51 this 
phase-shift contributes to the expression involved in the transversality conditions used to obtain diffusion. 


Remark 2.1. It is interesting to note that the invariant tori Fp^^py project onto the toroidal magnetic surfaces z = z* 
of the unperturbed ABC magnetic field in the configuration space T^. Moreover, the magnetic field on each surface 
is linear, i.e. ^abc\z=z* = (sin z*, cos z*, 0), and the trajectories follow the magnetic lines. Let us observe that the 
slope of the magnetic lines tanz* coincides with the ratio of the frequencies cf. Eqs. © and ®. 


2.2 Main Theorem: diffusion along a NHIM 

Let us consider the following Hamiltonian for the ABC system 

H = -^{Px — sinz — eC cos y)'^ + i^iPy ~ cos 2 : — eB sinx)^ -h -(pz — eC siny — eB cosx)^, (14) 

where we have introduced a scaling B = eB and C = eC. The following result states sufficient conditions for the 
existence of diffusing orbits: 

Theorem 2.2. Consider the Hamiltonian d of the ABC system with B > C D. Assume that the following 
hypotheses hold: 

Ai Considering the notation introduced in Section |Z 71 we define the functions M? = Mf{px,Py) as 

{py — COS z*) sin(x± -|- LOia) — (py — cos z^) sin(Fi -|- cuicr) — cos(Fi -|- wicr)^ da, 

:= C y ^{px -smz*)cos{y± + uJ2a) - {px - sin 2:°) cos(F2 -|-a;2cr) -p°sin(F2 -|-a;2cr)^ dcr, 
^(Py — cos z*) cos(x± -|- coia) — {py — cos z^) cos(Fi -|- ujia) + p^ sin(Fi -|- wicr)^ da, 

' J ^(Ps - sin2;°) sin(F2 -|- u:2a) - {px - sin z*) sin(y± -|- W2 (t) - p° cos(F2 -|- a; 2 < 7 )^ da, 


:= B 


MO 


:= B 


mO := C 


with Fi = Fi{a) and F 2 = F 2 {a), and where the notation x± (resp. y±) means that we take x_ (resp. y_) 
when we integrate in the interval (—oo,0), and x+ (resp. when we integrate in the interval (0,oo). We 
assume that there exists a non-empty set X = [oi, 5i] x {a 2 ,b 2 \, for positive values of a^, bi, such that and 
Mg do not vanish simultaneously, and the same for and M^, provided that (px^Py) £ F. 

A 2 Assume that for any value {px,Py) G X there exists a non-empty domain Jp^,py C with the property that 

X) := •Fp„.,py X {{px,Py)} C T X X 

iPx,Py)&I 
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is a domain, and when (x, y,Px,Py) S ^ there is a unique critical point t* = t*{x, y,Px,Py) of the map 
T 1 -^ All cos(x — ujit) + M 2 cos{y — 0 J 2 t) + M 3 sin(x — ujit) + M 4 sm{y — uj 2 t) , 
which defines a smooth function on D. 

A 3 Assume that we can chose a constant L > 0 such that for every (x, y,Px,Py) we have 

/ if\Px-Py\>L, 

\ A 1 A 4 - A 2 A 3 / 0, if\Pa;-Py\<L, 

where {Aj}j=i^ 2,3 cind {Aj}j=i^ 2 , 3,4 tire certain explicit functions depending on {x,y,Px,Py) that are defined 
in Section |3 cf Eqs. (I94l) - (I9^ and (I9^ - (I101I) . 

Then, given two pairs {Px,Py) £ i tind {Px,Py) G i and given 6 > 0, there exists e* = e*{d,I) such that if 
0 < |el < £* then there is a trajectory {x{t),y(t), z{t),px{t),py{t),pz{t)) of the system (fT4l) satisfying 

dist {{pI,pI), {px{0),Py{0))) < 6, 

dist {{pI,pI), {px{T),Py{T))) < 6. 


for some T > 0. 

We would like to emphasize that the above hypotheses are given in a very explicit way. To evaluate all the functions 
involved in the statement of Theorem 12.21 we only need to compute the coefficients {M?}j=i^ 2 , 3,4 in Hypothesis Ai, 
together with the partial derivatives t*, t*, T*y and Tyy, of the critical point in Hypothesis A 2 . As was sketched 
in the introduction, the proof of Theorem l2.2l consists in combining the internal dynamics on the NHIM with its outer 
(asymptotic) dynamics. Details are presented and discussed in Sections[3llU and|5] Finally, in Section| 6 ]we show that 
the hypotheses of the theorem can be rigorously checked in a computer assisted proof. 

Remark 2.3. Since the invariant tori 7'p^,py correspond to the toroidal magnetic surfaces of the unperturbed ABC 
magnetic field, c.f Remark I2.il Theorem 12.21 implies the existence of drift motions connecting any two magnetic 
surfaces (compatible with the set X) for the perturbed ABC system. This diffusion of charged particles is a very 
harmful phenomenon for the confinement of hot plasmas for fusion power generation, as explained in the introduction. 

3 Inner dynamics of the normally hyperbolic invariant manifold 

The study of normally hyperbolic invariant manifolds is a very classical (and important) topic and it has been ex¬ 
tensively considered in the literature. Most of the results that we use in this section are standard and can be found 
in ll2^ ITTII . Our purpose here is to present a basic overview, notation and perturbative formulas that we require to 
study the perturbation of the normally hyperbolic invariant manifold introduced in Section [2A] 


We recall that our goal is to study the Hamiltonian (fldl) for small values of e. Hence, we write H = 

He perturb- 

atively as follows 


He = Ho + eHi + £^H2, 

(16) 

where 


Ho = ^{Px- sinzf + ^{py- coszf + ^pl, 

(17) 


Hi = — C cos y{px — sin z) — B sin x{py — cos z) — Pz{C siny + B cos x), 

(18) 


^2 ^2 

H 2 = — + — + BC cos X sin y. 


( 19 ) 
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The unperturbed Hamiltonian Hq was studied in Section IZTl where we characterized the corresponding NHIM 
Aq. Now, we are interested in characterizing the perturbed invariant manifold together with the restricted dynamics 
on it (mainly the existence and approximation of invariant tori). To this end, we will follow closely the methodology 
introduced in the papers ifT^fTTlfT^ . 

Let us remark that the Hamiltonian (fT6l) is real-analytic. This will imply that all the objects obtained in this section 
will be of class C"", with arbitrarily large r (this follows from Fenichel rate conditions) so that we can omit all the 
discussions concerning regularity. This will simplify many technical issues, for example when applying averaging 
and KAM theory. The interested reader is referred to ifT^ [T8l for details on regularity. 


3.1 Normally hyperbolic invariant manifolds and perturbative setting 

Let M be a smooth finite dimensional manifold and let us consider a flow (jit, of class with r > 1, acting on M. 

Definition 3.1. Let A (Z M be a submanifold invariant under the flow, i.e., (pt{A) = A. We say that A is a normally 
hyperbolic invariant manifold (NHIM), if there exist a constant c > 0, expansion rates 0 < p < X, and a splitting for 
every x £ A 

T^M = © l ;" © T^A, ( 20 ) 

characterized as follows 


V & El 

\D(f)t{x)v\ < ce~^^ |u|, 

t > 0, 


V £ El ^ 

\D4>t{x)v\ < 

t < 0, 

(21) 

V £ T^A 

\Dft{x)v\ < u|, 

t G M. 



The classical theory of NHIMs guarantees that if A is normally hyperbolic, then it is persistent under small 
perturbations. Moreover, if the system depends smoothly on parameters, the manifolds —they may not be unique— 
can be chosen to depend smoothly on parameters. NHIMs are robust under perturbations, so we do not require a 
symplectic structure on M and ff Nevertheless, the problem considered in this paper is endowed with a symplectic 
structure and hence we will be interested in characterizing a symplectic structure on the perturbed NHIM. 

In order to apply the geometric mechanism for a priori unstable systems (c.f. we must compute explicitly 

some expansions in e of the NHIM associated to the Hamiltonian (fT^ . Notice that in our case we can model the 
NHIM by means of the canonical manifold x (see Section ITTI) . that is, we look for a parameterization 

: N ^ M, with Pe{N) = A^, characterized by the invariance equation 


X,oP, = DP,R, 


( 22 ) 


where Re is a vector field on N and Xe is the Hamiltonian vector field associated to He- Using the expansions 

Xe = Xo + eXi + + ..., 

Pe = Pq + sPl + £^P2 + ■ ■ ■ , 

Re = Rq + eRi + e^i?2 + • • ■, 


and equating terms in the expansion of s of the invariance equation (l22]) . we obtain (this approach was used in Uhl) 


0th order: Xqo Pq = DPqRq, 

1st order: (DXq o Po)Pi + Xi o Pq = DPq^i + DPii^o, 

2nd order: (DXq o Pq)P 2 + ^(D^Xq o Po)P^^ + (BXi o Po)Pi +X 20 PQ 

= BPqR2 + BPiRi + DP2^0) 

nth order: (DXo o Po)Pn - BPnRo - BPoRn = -Xn o Pq + Sn, 


(23) 

(24) 

(25) 


( 26 ) 
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where Sn is a polynomial in Xq, ..., Xn-i, their derivatives, Pq, ..., Pn-i, their derivatives, and Rq, ..., Rn-i- 
Clearly (see the discussion in Section IZTI) Eq. (|2^ has the solution 

Po{x,y,Px,Py) = {x,y,Z*,Px:,Py,pl) 

Ro{x,y,Px,Py) = iOi{px,Py)dx + U} 2 {Px,Py)dy 

where ui and 002 are given by (O and ([ 8 ]), respectively. In this case, since the unperturbed internal field Rq does not 
depend on the angular variables {x,y), the equations of the form ( 1 ^ lead to simple cohomological equations in a 
suitable frame. Hence, these equations can be solved explicitly using Fourier expansions. It is worth mentioning that 
there are more general theories that allow us to solve equations of the form (| 2 ^ even if the motion on the base is not 
quasi-periodic. 

As will be discussed in subsequent sections, the solution of equations (|2^ . (|2^ . (1^ . and (|2^ is not uniquely 
determined. We will use this freedom in order to obtain certain symplectic properties. More specifically, we follow the 
ideas in ifTTll to maintain the canonical symplectic structure on A^, so that we can easily characterize and manipulate 
the Hamiltonian associated to the restricted vector field 

3.2 Symplectic properties of NHIMs of Hamiltonian systems 

Let M be a symplectic manifold with symplectic form u, represented by a matrix-valued function ff, and let us 
assume that a C"" Hamiltonian Hq, with r > 2, has a NHIM Aq parameterized by Pq : M. Then, it is well 

known (c.f. |[2^lT7l l that for every perturbed Hamiltonian of class there exists a NHIM A^ parameterized by 
Ps of class (7^“^. Moreover, A^ is 0(e)-close to Aq in the sense. Here and in what follows, when we say 
that a map depending on parameters is of class C'' we shall mean that it is of class in all variables including the 
parameters. 

Given a family of Hamiltonians having a family of NHIMs A^ = P^{N), with P^ : N ^ M, we consider the 
maps Rg : N ^ TN corresponding to the vector fields restricted to the NHIMs. The maps P^ and R^ are related by 
the invariance equation ( 1221 ) . 

It is well known that the solutions of (l22l) are not uniquely defined, since we have the possibility of choosing 
different coordinates in the reference manifold N. It is natural to use this freedom to satisfy certain properties, like 
asking P^ to be a graph or asking R^ to be as simple as possible. In this paper, we are interested in choosing the 
solution that preserves the Hamiltonian structure of the problem, that is, we want that 

|:(P» = o. (27) 

The fact that this can be achieved was proved in IITtI . In this paper, since we need to perform some explicit computa¬ 
tions, we have to give some additional details on the procedure presented in ifTTl . The aim of this section is to explain 
the explicit computations required to handle a particular problem. 

A natural way to obtain (|27]) is to use deformation theory. Let us recall some standard definitions. Given two 
connected manifolds M and N, and given a family f^:N^M such that (x, e) 1 -^ feix) is in all its arguments, 
we define the infinitesimal deformation of as the vector field that satisfies 

^fe = PeOfe, 
de 

and we observe that = (^/e) o f~^ is defined on feiN) C M. 

Let be the infinitesimal deformation of the family P^ with initial condition Pq. It is clear that : A^ —)• TM, 
so we can consider the projections of according to the splitting (l 20 l) . Then we have the following result ifTTI : 

Proposition 3.2. Let us consider a family of parameterizations P^ : N ^ M with A^ = Ps{N). Assume that 
the infinitesimal deformation satisfies that the projection on the space T^jA^ vanishes for every x G A^. Then, the 
symplectic form is independent ofe, where is the original form u expressed in a basis of the splitting (1201) . 
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Proof. For the sake of completeness, we reproduce the proof given in ifTTl . First we observe that since u is invariant 
under the flow of H^, then also is and we have 




for every u,v £ T^M and t G M. Using the asymptotic properties in (l2TI) it is clear that u] = 0 if 

u G 0 £ and v G TxA,. (or vice versa). 

Then, using Cartan’s formula we obtain 


de 




where we used that is closed. Then, we have 






and we observe that the 1-form acting on u G T^N, is given by 


P*i'p^U^^e{x)[v\ = i'p^U^^e{Pe{x))[dPe{x)v\ = U>^^e{Pe{x))[Ve{Pe{x)),dPe{x)v\. 


By hypothesis, we have V^{P^{x)) G Pp^^x) e ® ^K{x) e dPe{x)v 

Peire^*,e{x) = 0 and we conclude that is independent of e. 


€ Tp^(.)A,. 


Hence, it must be 

□ 


Remark 3.3. A particularly interesting case arises if Aq = Po{N) is a NHIM for Xq and Eq. (1221) is solved per- 
turbatively. This is the situation considered in this paper. Property (1271) is important in order to have a canonical 
symplectic structure on A^, so that the averaging procedure (normal form) can be implement in the usual way. 


In the following we assume that M = N xT xW, with x M”, and we use the notation [u, Pu) £ X with 

u = (ui, ..., Un), Pu = {Pu,i, ■ ■ ■ ,Pu,n), and {v,py) G T X M. We endow M with the symplectic form 


which is represented by 




U 




n 

■ dpu,i A dui + dpv A dv, 
i=l 


(28) 


02nx2 


with 


Pr = 


On -Ir 


In 


Or. 


. 02 x 2 n 

where from now on we use the notation Onxm, Inxm, On = Onxn, and = Inxn, for the zero and identity matrices, 
respectively. Moreover, we denote by M^xn the space of m x n-matrices with real coefficients. 


Definition 3.4. Given a parameterization Pq : N ^ M of a NHIM, with W = T” x M"" and M = N x T x M., we 
say that Pq is compatible with the symplectic form u if 

I)Po{u,Pu)~^ Tln+lP>Po{u,Pu) = 

Similarly, we say that a frame 


<i-.Nxm?^+‘^ Tp^(N)M 

{u,Pu,i) -> {Po{u,Pu),CQ{u,Pu)i) ’ 

with Cq : N ^ M( 2 n+ 2 )x( 2 n-r 2 ). E symplectic if 


Oq{u,Pu) n-|-lOo( u,P n) — • 
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Let us also introduce some notation regarding derivatives that will be useful in computations. Given a vector field 
i? on a NHIM, and given a function ^ —>• M, we denote the Lie derivative of ^ with respect to R as follows 




2=1 


dn 




R^ 


'72+2 • 


(29) 


Moreover, given a parameterization P : N ^ M, and vector fields X and Ron M and N, respecfively, we infroduce 
fhe operafor 

7^p,x,i^(0 = DX o - Ljj(0 , (30) 

acting on funclions ^ —)• M. We exfend fhe nofafion in (l29l ) and (1^ componenf-wise for mafrix funcfions 

^ : N ^ Mmxn- In ofher fo simplify fhe nofafion, we will wrife TZq = TZpq^Xo,Ro- 

Given a parameferizafion Pq : N ^ M of a NHIM, wifh N = x and M = A^ x T x M, we can fake 
derivatives af bofh sides of fhe invariance equafion ^ 00^0 = DPq^o thus obfaining 


DXo o PoDi^o = D(DPoPo) = Li{g(DPo) + DPqDPo- 
This means fhaf fhe fangenf vectors of Po(At) partially characferize fhe action of fhe operafor TZq in (l30l) as 


'^o(DTb) = DPqDPo- 


Since Po{N) is normally hyperbolic, fhere exisf maps Wq : N —)• M(2n+2)x2 parameterizing fhe normal bundle of 
Po{N), and Lq : A^ —)• M 2 X 2 such fhaf 

RoiWo) = WoTo. 

From now on, we assume fhaf Fq is diagonal, and due fo fhe Hamiltonian strucfure we can wrife 


Moreover, if we assume fhaf Pq is compatible wifh fhe symplecfic form u, fhen if furns ouf fhaf fhe mafrix Wq can be 
scaled in such a way fhaf fhe juxfaposed mafrix Cq := (DPq Wq) G M(2n+2)x(2n+2) defines a symplecfic frame as in 
Definition 13.41 

The operafor TZq infroduced above appears in fhe perfurbafive equations (|2^ - (|2^ obfained in Secfion 13.11 The 
following lemma approaches fhe sfudy of fhese equafions using fhe previously consfrucfed frame. If is worth men- 
fioning fhaf fhe facf fhe frame C is assumed fo be symplecfic is nof really necessary. Nevertheless, if simplifies some 
compufafions (for example fhe compufafion of fhe inverse Cq^). 

Lemma 3.5. Assume that Pq : N ^ M satisfies Xq o Pq = DPqPq, with A^ = T" x M” and M = A^ x T x M. 
Given a map rj : N ^ we consider the following equation 


DXq o Pq^ - D^Po - DPqp = r] (31) 

for the unknowns ^ : N ^ and p : N ^ Then, using the symplectic frame Gi associated to the matrix 

Cq = (DPq Wq) constructed above, it turns out that Eq. (1311) leads to 

-Lijo (1^) + DPol^ = f + p (32) 

-Lfi,(p) + rop = 7}^ (33) 


where ^ 

i = Gol = DPqI^ + WQi^ and fj = = -Ttn+iC^^ Dn+iti , 

with fj^ : N ^ + : At ^ ■, N ^ and : N ^ M^. 
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Proof. Let us observe that the fact that C is chosen to be symplectic allows us to compute the inverse of Cq as follows 

Cq — rifi-i-i. 

We also notice that the action of TZq on the matrix Cq^ takes the form 

'^o(C'ol) = TZo{Co)i - CoLijpd), 

and that 


Cf^noiCo) = - Qn+lC^Qn+l (DXo O PoCo - Lr,{Co)) 

^DPjQn+i^PoBRo DPo^fl„+ilLoro\ 

WjQn+iBPoBRo Wj 0„+ikLoro J 


— ^n+l 


Di?o 02nx2 
02x2n To 


Di?o 02nx2 

02x2n To 


Introducing ^ = Cq^ = into Eq. (I3l]l, we obtain 


Di?o 02nx2 

02x2n To 


|C 


L/Jo 


- Cf^PPoP = 


Then, we observe that 


-Cf^PPoP = ^n+lC^^n+l^PoP = ^n+1 [wjnn+i 


fP>Po fln+lDPo\ 
DPo y 


p = 


hn 

02x2n 


Finally, using the symplectic structure, we introduce and p^ as in the statement of the lemma, thus ending up with 
the equations (l32l) and (l3^ . □ 

It is standard to check that the solution of Eq. (l3^ is unique. In our particular case (see computations in Sec¬ 
tion (33]), it turns out that Rq produces an integrable quasi-periodic motion in N, and hence, we can solve (1331) using 
Fourier series. In particular, if we have a function /3 : —)• M expressed in Fourier series as 


/3{u,Pu) = cos{k • u) + iSf^ipu) sm{k • u)) , 

fceZ" 

with /3 q“ = 0 , then it turns out that the solution ^ of the equation — L/jp = /3 is given by 

C{u,Pu) = i^k°''(.Pu) cos{k • u) + ?r(p„) sin(A; • u)) , 


with 




+ UJ-& 


Sin 

k 


ef = 


prx - w 


ocos 


A^ + (cu • k)'^ ’ A^ + (cj • k)'^ 

In case that Rq takes a more general form, Eq. (l3^ can be solved using the asymptotic properties of the cocycle. 
As was mentioned in Section [XT] the solution of Eq. (I32l ) is not unique. A simple choice consists in taking 


(34) 


i^ = 02nxl, p=-p^, 


(35) 


but, in general, this solution will not determine a parameterization which is compatible with the symplectic structure 
of the problem. The final goal of this section is to compute the deformation of the symplectic frame ^ with respect to 
the perturbation pai'ameter and to combine Proposition 13. 2l and Eemma lXSl in order to obtain the canonical symplectic 
structure in the deformed NHIM. 
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Assume that : N ^ M, with = T” x and M = A^xTxM, isa family of parameterizations satisfying 
X^o where is a family of Hamiltonian vector fields with the symplectic form u: given by (1^ . Let 

us consider Pg, the infinitesimal deformation of the family Pg with initial condition Pq. A simple computation shows 
that 

7/7-j 

— 7 ^ = Pq + 2P2e + SPse^ + ... = Pq o Pq + (DPq o PqPi + Pi o Po)e + ... , 
ae 

thus obtaining 

0 th order: Pi = Po o ^0 

1st order: 2 P 2 = Pi o Pg + DPg o PgPi 

nth order: nP„ = Pn ° Po + Sn 

where Sn is an explicit expression depending recursively on the previously computed objects. 

Let us consider the first order correction determined by Eq. (I24b . We apply Lemma 1X5] with 

^ = Pi, P=Ri, r] = -XioPQ 

and we consider the unique solution of Eqs. (1^ and (l3X) satisfying Eq. (1^ . In Eq. (l3^ we observe that Pi is 
proportional to Pg. Hence, it turns out that the deformation Pg vanishes on the central directions. By Proposition 13.2[ 
we conclude that the reduced vector field Pi is a Hamiltonian vector field wifh respecf to fhe form 17°. 

The second order correction is nol so simple. On fhe one hand, we observe fhaf P 2 and Pi are no longer propor- 
fional. On fhe ofher hand, we have to consider Proposition 13. 2l on fhe deformed symplecfic frame. Eet us assume fhaf 
we have compufed Pi, Pi, and also fhe first order correction of the symplectic frame, that is, Cg = Cg + eCi + O (e^). 
Then, we express the infinitesimal deformation Pg on the frame Cg perturbatively as 

C'-'Pg(Pg) = 

C'g-^Po O Pg + e(Cg-^DPg O PgPi + Cg'^Pl O Pg - Cq^CiC^^Po O Pg) + 0(^2). 

By construction, it is clear that 

Cg-lPo O Pg = Cg-'Pl = Cg-^Coll = |l = ') ■ 

We ask the same condition for the e-order terms, thus obtaining that 

Cg-^DPo O PgPi + Cg-lpl O Pg - Cg-^ClCg-lpg O Pg = ^ 

for certain (^ : iV —)• whose expression is not important for us. Then, we use again that Pg o Pg = Pi = Cq^i, we 
replace Pi o Pg using (iTTl) . and we write P 2 = Co^ 2 , thus obtaining the condition 

26 - Cg-^Ci ') = (39) 

that determines the first 2n components ^2 of ^ 2 - Therefore, we can solve the second order correction of the invariance 
equation, given by (|25]) . using Eemma lXSl with 

C = P 2 , P = R2, 7? = -X 2 O Pg + DPiPi - ^D^Xg O PgPf® - PXl O Pi 

and choosing the unique solution obtained by fixing ^2 satisfying Eq. (l39l) . Then, fhe corresponding correction of fhe 
reduced vecfor field, 

P 2 = DPg^^^ - L,J0 (e^) - fi^, (40) 


(36) 

(37) 

(38) 
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is a Hamiltonian vector field with respect to the form . 

Finally, we need to give a simple recipe to compute the first order correction Ci of the symplectic frame. The 
construction is analogous up to any order, but this is enough for our purposes. We will construct the frame taking 
Cl = (DPi Wi), where Wi is computed as follows. On the one hand, we assume that we have computed = 
Pq + ePi + C(e^) so that we have (the computation is direct) 

'Ppe,Xe,Re{^Po + = (DPq + eDPi)(Di?o + eDPi) + O(e^), 

where we recall that Tlp^^Xe,Re given by Eq. (l30l) . On the other hand, we look for Wi and Fi is such a way that the 
action of Tlp^,Xe,Re ih® matrix Wq + eWi is given by 

= {Wo + eWi){ro + eFi) + ©(e^). 

We observe that this condition is satisfied if 

(DXo o Po)Wi - - WoW - WiTo = -Si, (41) 

where 

5i := Lr, {Wo) - BXi o PoWo - B^Xo o PoPi ® IFo ■ 

Again, fhe solutions of Eq. (IdTl) are obfained by considering fhe action of fhe unperfurbed operator TZo- In the 
following resulf, analogous fo Eemma lTSl we study the above equation. 

Lemma 3.6. Assume that Pq : N ^ M satisfies Xo o Po = BPoRq, with At = T” x M"" and M = At x T x M. 
Assume that the pair Pi and Ri is a solution of equation (1241) . that is, we have 

{Xo + eXi) o {Po + ePi) = (DPq + eDPi)(Po + eRi) + 0{e^). 

Then, using the symplectic frame (£ associated to the matrix Co = (DPq Wo), it turns out that equation (1411) leads to 

-Lr, {W^) + BRoW^ - W^To = , (42) 

-Lr, (Htf) + ToWi^ - Wi^Fo = -Sf - El, (43) 

where 

Wi = CoHti = DPoVFf + WoW^ and = -On+iCj" B^+iSi . (44) 

Proof We recall that the frame C satisfies 

no{Co) = DXo O PoCo - B{Co)Ro = Co {. 

\tX2x2n to / 

Then, we compute the action of Pq on Wi = CqWi as follows 

no{CoWi) = Co Wi - CohRfiWi), 

\B2x2n to / 

and we introduce this expression into (|4TI) . thus obtaining 

Co ( #1 - CoLr,{Wi) - kFori - CokFiEo = 5i. 

\^2x2n to / 

Using the symplectic properties of the frame, we multiply both sides by Cf^ = —U„+iCj^fin+i and we end up with 
( 1^1 - D(lUi)Po + fln+iCo^fln+iH^ori - IE,To = -U„+lCo^ Bn+lSl . 

\^2x2n to / 
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Finally, we observe that 


{wjnT^w:) = {°T) ' 

and using the notation in (l44l) we obtain the equations (|4^ and (|43]) . □ 

Finally, we discuss the solution of equations (|4^ and (|4^ . On the one hand, we observe that equation (|4^ is 
similar to equation (1^ in the sense that it can be solved using Fourier series, obtaining a unique solution. On the 
other hand, we observe that the diagonal part of the left hand side of equation (|4^ is resonant. We can avoid this 
resonance by selecting Fi. To this end, we consider the particular choice 

Fi = diag(5f), 

where (•) stands for the average with respect to the variables u £ T". Obviously, this choice preserves the diagonal 
character of the matrix Fg = Fq + eFi + .... 


3.3 Perturbative computation of the NHIM of the ABC system 

The goal of this section is to compute the NHIM Ag associated to the ABC system in the perturbative setting given by 
Eqs. (fT^ - ([T^ . We follow the notation and methodology described in Section [TH 

First, it is convenient to reorder the phase-space coordinates as (x, y,Px,Py, z,pz) rather than (x, y, z,px,Py,Pz)- 
In analogy with Section we have {u,pu) = {x,y,Px,Py) £ iv = X and = iz,pz) £ T x M. Then, 

we consider the symplectic form u, and its matrix representation given by 


u = dpx A dx + dpy /\dy + dpz A dz, 




/ 04x2^ 

V02X4 Jd? ) • 


With the above notation, we have 

We start with the explicit characterization of the unperturbed problem, giving rise to the expressions 


/ Px- sin z 
Py — COS 2 


Pz 

\Px cos Z-Py 

and the corresponding derivatives 


\ 


/ 

X 

\ 




y 



, ^^0 = 


Px 




Py 





= arctaia{px/py) 

+ TT 

z) 


1 

Pl = 0 

/ 


/O 

0 

1 

0 

— cos z* 

o\ 


/I 

0 

0 

0 

\ 

0 

0 

0 

1 

sin z* 

0 


0 

1 

0 

0 


0 

0 

0 

0 

0 

0 

, DPo = 

0 

0 

1 

0 


0 

0 

0 

0 

0 

0 

0 

0 

0 

1 


0 

0 

0 

0 

0 

1 


0 

0 

Py\-^ 

-PxX~ 

-4 

Vo 

0 

cosz* 

— sin z* 

A2 

0/ 


\0 

0 

0 

0 

/ 


where we recall that A = sinz* = —pxjX'^, and cosz* = —py/\^ (see computations in Section l2n) 

Notice that the parameterization Pq given above is compatible with the form u> according to Definition 13.41 

To obtain the unperturbed symplectic frame we take the columns of DPq and we complement them with the 
eigenvectors of DAq o Pq of eigenvalues A and —A, that we suitable scale in order to obtain a symplectic frame. 
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Specifically, we take 


C-o 


no 0 
0 1 0 

0 0 1 

0 0 0 

0 0 

\0 0 0 


0 {V2l2)pyX-y^ 

0 -(V2/2 )p^A-^/2 

0 0 

1 0 

-p,A-4 (V 2 / 2 )A-V 2 

0 (V 2 / 2 )A 1/2 


{V2l2)pyX-y^ \ 
-(V2/2)p^A-7/2 
0 
0 

-(V 2 / 2 )A -^/2 

(V 2 / 2 )Ai /2 


and we left as an exercise to the reader to check that C'o(x)^fl 3 C'o(x) = fls, where tfa is the matrix of the canonical 
symplectic form. The inverse of Cq is given by 


Cn-l = 


/I 

0 

0 

0 

0 

Vo 


0 

0 

1 

0 


0 

0 

0 

1 


0 

0 

0 

0 


-PyX-^ 
PxA ^ 
0 
0 


-(V2/2)pyA-7/2 (^/2 )p^A-^/2 (V2/2)A1/2 {V2l2)X-^/^ 

{V2l2)pyX-y^ -(V2/2)p^A-7/2 -(V2/2)AV2 {^j2)X-^/^J 


and it turns out that this frame allows us to reduce DXq o Pq as follows 


Cq^BXo o PoCo = 


/O 0 pIx-^ + 1 

0 0 -PxPyX~^ 

0 0 0 

0 0 0 

0 0 0 

Vo 0 0 


-PxPy>^~^ 0 0 \ 

pIX-^ + 10 0 

0 0 0 

0 0 0 

0 AO 

0 0 -A/ 


Notice that this expression corresponds to Cq^TZo{Co) since in this case C'(^^D(C'o)Po = 0. Finally, the reduced 
vector field is given by Rq = ojidx + uj 2 dy, where we recall that wi = px{l + X~^), and U 2 = Py{l + A“^). 

To obtain the corrections Pi and Pi of the parameterization and the reduced vector field, respectively, we consider 
the equation 

(DXo o Po)-Pi — D-Po-Ri — = V 

where 

C cos y \ 

B sinx 
— BuJ 2 cos X 
CiOl siny 
C sin y + B cos x 
cos z* cos y — B sin z* sin x) 

Using Lemma [+31 with Pi = CqC and Pi = p, we obtain the equivalent system of equations 

(PyA"® + 1)|3 - PxPyX~^U - LiJo(li) = 7}i + Pi , 

—PxPy^ + (Px^ ^ + 1)|4 — Lijp(|2) = t72 + P2 ) 

-L/Jods) =fj3+P3, 

-^RoiU) = m + Pi, 

AIs - L/Jods) = 7)5 , 

-Ale - Lijpde) = % > 


(45) 

(46) 

(47) 

(48) 

(49) 

(50) 


V = -Xi o Po = 


V^ 
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where = uAii + uj 2 dyii, and 


/m\ 


/ Ai cos y + A2 sin x \ 

m 


A3 cos y + A4 sin x 

m 


A5 cos X 

m 


Ag sin y 

m 


A7 cos X + Ag cos y + Ag sin x + Aig sin y 

KmJ 


\Aii cos X + Ai 2 cos y + A13 sin x + A14 sin y/ 


The coefficients Aj, i = 1,... 14, are functions depending on the action variables Px,Py, given by 


Ai — (7(1 + PyX 

A2 = - BpxPyX~^ 

^3 = - CpxPyX~^ 

a, = b{i+pIx-^) 

A5 = — BuJ 2 

Ag = Cloi 

Ay = ^/2(A® + A2p2 +p 2 ) 5 ^-ll /2 


As = - V2l2CpyX-^^^ 

Ag = V2l2BpxX-^/‘^ 

Aio = V2/2{X^ + A2p2 +p2)^^-ll/2 

An = - v/2/2(A® + XVy + pl)BX-^^/^ 
Ai2 = - V2l2CpyX-^/^ 

Ai3 = V2j2BpxX-^^^ 

Ai 4 = - V2/2(A® + A2p2 +p2)(J^-ll/2 


(51) 


The solution of Eqs. (|4^ and (l50l) . using Fourier series, is obtained using (IMl) 

^5 = ill cos a; + il 2 cos y + B^ sin x + ^4 sin y , 
^6 = ill cos X — B 2 cos y — B^ sin x + il 4 sin y , 


where the coefficients Bi have the following expressions: 


Bi 


A7A + W1A9 _AsX + U}2Aio _ AgA — W1A7 _ AloA — W2^8 
X^ + ul ’ A2+a;i ’ X^ + cof ’ A^ + a;2 


which are functions depending on the action variables Px,Py The solution of Eqs. (1431)- (1431 ) is given by ^1 
I 3 = I 4 = 0 and 


ill = p = 


/—Ai cos y — A2 sin x\ 
—A3 cos y — A4 sin x 
—A5 cos X 

\ -Agsiny / 


6 = 


(52) 


By construction, the vector field ili in (l52l) is Hamiltonian with respect to the symplectic form dpx A dx + dpy A dy 
(c.f. Section lT2l) . 

The specific computations regarding Ci, R 2 and P 2 are omitted, since they will not be used in what follows. 
The only thing that we need to know in the next section is which resonances appear in the averaging process of the 
Hamiltonian corresponding to R 2 . We remark that in our problem, it turns out that R 2 is a trigonometric polynomial 
of degree 2. This claim follows from the the construction explained in Section 13.21 and the fact that we know the 
degrees of the functions Xq, Xy, X 2 , Po, Pi and Cq that appear in Eqs. (l39l) and (l40l) . 


3.4 Invariant tori on the NHIM 

From the computations presented in Section 1331 we obtain that the dynamics reduced to the perturbed NHIM is given 
by the Hamiltonian system: 

re{x,y,Px,Py) = ro{px,Py) +ri{x,y,px,Py)e + r 2 {x,y,px,Py)e‘^ + 0{s^). 


(53) 
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The Hamiltonian functions satisfy Ri = —where Ri is the reduced vector field on the NHIM computed in 
Section |33] Specifically, we have 

ro{Px,Py) = 

ri{x,y,Px,Py) = Asinx - Aecosy , 

r 2 {x,y,Px,Py) = -4i5 + Ai6 COS X + An COS y + Aigsinx + Aig siny 

+ A20 cos(2x) + A21 cos(2y) + A22 cos(x + y) + A23 cos(x - y) 

+ A24 sin(2x) + A25 sin(2y) + A26 sin(x + y) + A27 sin(x - y), 

where A 5 and Aq are given in Eqs. (ISTT) and fhe remaining coefficienfs are certain explicit functions of {px,Py) whose 
explicit expressions are not important in the computations performed later. Here we are denoting as {x,y,px,Py), 
with abuse of notation, the reduced variables on the perturbed NHIM, but they are not the same as the coordinate 
variables in the phase space x However, at first order in e the parameterization is a graph (see Eq. (l36l) l. and 
hence, the reduced variables and the coordinate variables only differ in terms of order 


3.4.1 The global averaging method 

The task now is to characterize invariant tori on the perturbed NHIM. The idea introduced in ifT^ consists in perform¬ 
ing several steps of averaging in a global way on the whole NHIM. To this end, a normal form procedure is applied but, 
when we are close to a given resonance, the resonant normal form is defined by evaluating fhe corresponding coeffi- 
cienf on fhe resonanf manifold (see also HU). If is worfh menfioning fhaf since fhe problems considered in 
are non-aufonomous, a suifable projection is fhe so-called projection along fhe fe-direcfion. In our case, due fo fhe facf 
fhaf fhe sfudied Hamiltonian is autonomous, the orthogonal projection is more appropriate to perform computations. 

Although we are interested in the ABC system, the discussion of this section is presented in a general setting. 
This allows us to use a more convenient notation and, moreover, we think that it will help the reader to link with the 
exposition in ifT^ HU and to consult the details that we omit in our discussion. 

Eet us consider a Hamiltonian system on x of the form 

m 

h{u,pu) = ho{pu) + '^e'-hi{u,pu) + , (54) 

where every hi is written in Eourier series as 

hi{u,pu)= ^ {hl°^{pu)cos{k-u) + hf^{pu)An{k-u)) , (55) 

where Zi is the support of the Eourier series, which is assumed to be a finite set. Eor the sake of consistency we take 

Lsin _ n 
''' 1,0 — 

The averaging of Eq. (l54l) consists in performing (recursively) a suitable change of variables in such a way that we 
obtain a new Hamiltonian system depending on the variables rx G in a simple way. Setting {h)o{u,pu) := h{u,pu), 
let us assume that we have performed m — 1 > 0 steps of averaging, so we have 

m—1 

{h)m-i{u,Pu) = ho{pu) + ^ e%{u,pu) +e"^hm{u,Pu) -hC>(e™+^). 

i=l 

Then, given a Hamiltonian system e^gm with time-1 flow we infroduce fhe new Hamiltonian 


{h^raiAi Pu) 


O (lA'^{u,Pu) 
m—1 

ho{Pu) + ^ e^hi{u,pu) -h e”" 


hm Pu) + {ho,gm}iu,Pu) 


+ 0 (e™+i), 
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and we ask it to be as simple as possible by taking 

hm{u,Pu) + {ho,gm}iu,Pu) = hm{u,Pu) ■ 
Here {•, •} is the Poisson bracket, defined as 


Sf 1 ( 9f dg df dg \ 

^^\^'^idPu,i dPu,idui) ' 


Using an expansion as in (l55l) . we obtain the following set of equations for the Fourier coefficients: 



(56) 


for every k G Z”\{0}, and we take g = g, and /i^“g = 0, so that p^g and gr^“g can take any value. In these 
equations w = a;(p„) 

Definition 3.7, Given a Hamiltonian h : I G M” —)■ M,/or each k G Z’^\{0} we define the resonant set 


Rfc := {Pn G Z : uj{pu) ■k = 0}. 


Let us assume in what follows that the function ut {pu ) • k has no critical points on , so that the resonant manifolds 
are smooth surfaces (a condition that depends only on the unperturbed problem and that is certainly satisfied by the 
ABC system). Then, it makes perfect sense to introduce some additional definitions and notation. Indeed, given 
a resonant set Rfc and a small constant L > 0, we denote the tubular neighborhood of R^ of radius L (measured 
with the Euclidean metric) as Tub(Rfc, L). Moreover, for every resonant set we introduce the orthogonal projection 
Hfc : Tub(Rfc, L) C —)• Rfc. Finally, given a resonant set R^ we denote by dist (p„, R^) the Euclidean distance of 
the point pu to the manifold R^. 

Notice that R^ = Krnk for any m G Z. Then, given two sets R^ and R^, we have, generically, the following 
trichotomy: 

• They are the same manifold: R^ = R^, i.e, k = mi for some m G Z. 

• They do not intersect: R^ n R^ = 0. 

• They intersect transversely in a manifold of codimension two without boundary. 

It is worth mentioning that the third case does not play an important role in our problem. Indeed, for the ABC 
system, resonant sets are 1-dimensional manifolds and their intersections define sets of dimension zero. The case 
of higher dimensions has been discussed recently in ifT^ proving that the existence of multiple resonances is not a 
limitation to prove diffusion. 

If there is a finite number of resonant sets, it is clear that we can choose a constant L > 0 small enough such that for 
every pair A;, £ G Z” we have either R^ = R^ or R^ n Tub(R£, L) = 0. Then, following ifT^ . we construct a solution 
of Eq. (l56l) in a global way, that is, for all values pu G X. Of course, we only want to modify the Fourier coefficients in 
the support of the series, that is, ifk^ we take = h^^^kiPu) = 0 , and hence = 0 - 

Then, if A: G Zm, we take 



where : M — )■ M is a fixed C°° function such that = 1, if t G [—1, 1], and V'(A) = 0, if A ^ [—2, 2]. The Fourier 
coefficients of the Hamiltonian g^. are obtained from Eq. (l56l) . passing to the limit when p„ tends to Rfc. For details, 
we refer to Lemma 8.8 in ifT^ and to Lemma 10 in ifT^ . With this choice we distinguish two different zones: 
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• Non-resonant region: If ^ Tub(Rfc, 2L), we have = 0 = h^^^^{pu) = 0. 

• Resonant region: Ifp„ e Tub(Rfc,L), we have h'^^%{pu) = h'"m%0^k{pu)), and Tf^^kiPu) = K^,k(^k{Pu))- 

Remark 3.8. The choice of L is arbitrary. This implies that we do not need to study the regions at a distance between 
L and 2L of the resonant set R^. 

3.4.2 Adapted coordinates on the averaged system 

Let us now apply the averaging procedure described in Section [3.4. II to the reduced Hamiltonian (l5^ . In this case, 
resonant sets are expressed as 

Rfc = {{Px,Py) € T : Wifci + 0J2k2 = 0} 

where k = {ki,k 2 ), the set Z C {px > 0} x {py > 0}, and the frequency uj = (cji, a; 2 ) is given by Eqs. (|7]l and dUl. 
Then, it is clear that there are no resonances associated to the averaging of order |A:| < 1, since Px / 0 and py / 0, 
and so we have cji / 0 and UJ 2 / 0. For the same reason, in the averaging of order |A:| < 2, we must deal only 
with the set tui — a ;2 = 0, that corresponds to the straight line px = Py The orthogonal projection associated to this 
particular resonance, that we simply write as H, has the following explicit expression: 


n(P.,P.) = (^.^). (57) 

• Non-resonant region: we can eliminate all the terms in ri{u,pu) and r 2 {u,pu), so the second order averaged 
system is given by 

{re) 2 {x,y,Px,Py) = ro{px,Py) + 0{e^). 

Neglecting the 0{e^) terms, we obtain an integrable unperturbed system. The invariant tori of this unperturbed 
system are given by the level sets 


Py = £2- 

When the perturbation 0{s^) is taken into account, KAM theorem guarantees that most of these invariant tori 
persist for the perturbed system, covering the non-resonant region except for a set of measure of order 
We remark again that, since our problem is real analytic, we do not need to care about the technical difficulties 
regarding regularity in the KAM theorem. 

• Resonant region: we consider the projection (1571) . and we obtain that the second order averaged reduced 
Hamiltonian is given by 


{re) 2 {x,y,Px,Py) = ro{px,Py) +£ A 23 (n(p^,py)) cos(x - y) + A 27 {U{px,Py))sm{x - y) +0{£^). 


Then, it is natural to perform a canonical change of variables in order to introduce a resonant angle: 


61 =X, h=Px+Py, 

d2 = x-y, I 2 = -Py, 


thus obtaining the Hamiltonian 

{r£)2i0i,02,h,l2) = roih + I 2 , -h) + £^(^2z{^, ^)cos(6»2) + ^ 27(^7 + 0{£^). (59) 
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The next step is to perform a Taylor expansion around the resonance. It is clear that the resonance px = Py is 
equivalent to I 2 = —1\I2. Hence, we consider the expansion I 2 = —/i/2 + b and we write the unperturbed 
Hamiltonian as follows 


ro(^ + 4 - 5) _ ro(4, 4) + i (1^ - 2^^ + + 0(^3)^ 

where we have used that wi = a ;2 at 5 = 0. Moreover, using the specific expression of tq, it turns out that we 
can write the Hamiltonian (1591) as 


(?’£) 2 ( 6 'i, 6 ' 2 ,/i,-y+( 5 ) = ^ + ^+(l+V2 /^ t)cos(6»2)+A27(t, 

modulo terms of order This corresponds to a pendulum-like Hamiltonian system in the variables 

(02, (^) depending on the variable l\. In other words, we observe that Ii is an integral of motion of the truncated 
Hamiltonian (l60l) and the motion of the variables (02, b) is described by the system 


02 = 2(1 + V2/f ^)<5, 

b = £^ ^A23(yi sm( 02 ) - ^27(t) cos( 02 )^ 


The above system has a hyperbolic equilibrium point at (02, b) = (02(/i), 0), and we denote by H* the energy 
of this point. Then, the level sets of Ii and {ri .)2 are characterized as follows 


.^1 = ei, 

+ (1 + \/2/^ ^) b ‘^ + ^A23(-y, y) cos(02) + A27('^, ^;) sin(02)^ + 0 {£^, (5^) = 62 . 

We observe that these level sets have different topology depending if 62 > H* (two primary tori), if 62 = H* 
(two whiskered tori with coincident whiskers) or if 62 < H* (two secondary tori). Again, applying the KAM 
theorem to consider the effect of the perturbation, we obtain that many of the invariant tori in the previous 
picture persist, covering the resonant region except for a set of measure Cl(e3/2) We refer to llT^lT^lT^ for 
full details on the application of the KAM theorem close to the separatrix. As before, we do not discuss here 
the specific technical details since they are covered by the fact that our Hamiltonian is real-analytic. 


We have obtained an approximation of the level sets that characterize the invariant objects inside the NHIM. Such 
level sets are not written in terms of the original variables of the problem but in the averaged variables. In the following 
result we translate the previous construction into the coordinate variables in phase space. 

Proposition 3.9. Let us consider the original Hamiltonian system 

H, = Ho + eHi + e‘^H2, 

where Hq, Hi and H 2 are given by Eqs. (QT])- ([Ill). Then, the invariant tori inside the NHIM are characterized by the 
level sets 

Px - eB^sinx + 0{e^) = ei, 

Py - y + , 

in the non-resonant region, and 

Px+Py -eB^sinx - e(7^ cos y-h O(e^) = ei, 

+ (1 + V2e^ ~ Py cosyj C>(e^) = 62 , 

in the resonant region. Recall that the frequencies uji and UJ 2 are defined in Eqs. (jT]) and (H)). 


(63) 
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Proof. We only have to undo the different changes of variables in the averaging construction previously explained. In 
particular we recall that we defined 

{rf)i{u,pu) = reO cj)3i(^u,pu), 

where is the time-1 flow of a Hamiltonian egi satisfying ri -|- {ro, } = 0, with ri = A 5 sin x — Aq cos y. The 

expressions of A 5 and Ag are given in Eqs. (BTI) . Since there are no resonances involved, we can solve the previous 
equation taking 

gi {u, Pu) = Gi cos X + G 2 sin y, 

with Gi = and G 2 = —C^. This means that we have to invert the change of variables 

xr^ x + edp^gi + O(e^), 
y^y + edp^gi+Ois"^), 

Px ^Px -edxgi +0{e ), 

Py^Py- ^dygi + 0{e^), 

that lead to the averaged system. As explained at the beginning of Section lT4l the reduced variables are the same as 
the coordinate variables up to terms of order so we can safely assume that the variables {x,y,Px,Py) in Eq. (I64l ) 
are the phase space coordinates. 

Eet us first consider the non-resonant region. The unperturbed invariant tori of the averaged system are given by 
the level sets px = ei and Py = 62 , cf. Eq. (15^ . Inverting the change of variables (IMl) we obtain that the surviving 
invariant tori satisfy the expression in Eq. (I62l ). 

In the resonant region, we obtained that the unperturbed invariant tori of the averaged system are given by the 
level sets in Eq. (IMT) . Eollowing ifTSl . we first replace the expression Ii = ei into the second expression in (IMT) . thus 
obtaining the equivalent system 


h = ei, 

| + ^ + (l + V2er‘)i^ + 0(E^) = e,. 

This choice will simplify subsequent computations. Then, we recall the definition of the variable 5 and we invert the 
change of variables {px,Py) thus obtaining 

5 -I 2 + - --Py+-. 

Then, inverting the change of variables (1^ . we obtain Eq. (1^ . □ 


We would like to remark that the terms of order e in Eqs. (1^ and (l6^ will be important in the computations of 
Section [21 These terms are not necessary in due to the fact that the unperturbed outer dynamics is the 

identity and hence there is no phase-shift. 


4 Outer dynamics of the NHIM 

In this section we consider the outer dynamics of the NHIM for the perturbed system. This dynamics is modelled 
by the so-called scattering map of a normally hyperbolic invariant manifold with intersecting stable and unstable 
invariant sets along a homoclinic manifold. This remarkable tool was introduced in ifTSl to study Arnold diffusion in 
the context of periodic perturbations of geodesic flows in T^, and it was crucial for applications in |[T^[T 8 ll 20 ll 26 l . 
The paper ifTTl contains a complete description of the geometric properties of the scattering map, together with a 
systematic development of perturbative formulas for its computation. 

In Section l4fil we recall the construction of the so-called Melnikov potential, which was introduced in ifTSl . in the 
setting considered in this paper. In Section we present a brief definition of the scattering map and we obtain its 
first order approximation for the case of the ABC system. 
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4.1 The Poincare-Melnikov function 

As was discussed in Section[3l for small values of s there exists a perturbed NHIM, denoted by A^, together with local 
invariant manifolds WjQj,(Ae) and VFi“^(Ae). These manifolds are 0(e)-close to A and 14^®(A) = 14^“(A), respect¬ 
ively. As usual, we globalize the invariant manifolds as 44^®(Ae) = Uf<o '/’f (^ioc(^£))’ = Ut>o (^ioc(^£))’ 

where cpf is the flow of the perturbed Hamiltonian The intersections of the stable and unstable invariant manifolds 
are given by the following proposition. All along this section we use the notation introduced in Section l2dl for the 
unperturbed problem. 

Proposition 4.1. Let us consider an analytic Hamiltonian system of the form H^{q,p) = HQ{q,p) + £h{q,p,e), 
having a NHIM A^, where the unperturbed Hamiltonian Hq is given by Q. The homoclinic intersections of the 
invariant manifolds VF®(Ae) and 14^'^(Ae) are described, at first order in e, by the critical points of the Poincare 
function (also known as Melnikov potentialj.- 


/ OD 

- Kfliu* +u±),0)da, 

■OO 


( 66 ) 


where is the time-a flow of the unperturbed Hamiltonian Hq. In particular, {u* + u ±) is given by Eqs. (1101) 
and (1131) . and is given by Eq. dill) . Recall that the compact notation u± means that we take Uj^for cr > 0 and 

U- for (T < 0. 

We observe that this expression of L{T,x,y,px,Py) differs from the one used in ifT^ [T^ [T^ by the fact that 
it depends on the phase-shift. A Melnikov potential of this type is given in Proposition 3 of BTIl and analogous 
expressions can be found in ifldlfTSl . We invite the reader to compare Proposition 14. 1 1 with Theorem 32 in ifTTl that is 
stated in a more general setting. For the sake of completeness, we present here a complete proof of this proposition that 
may be of valuable help for the general reader. The arguments, which we adapt from ifT^ . are standard in Melnikov 
theory and well known to experts. 

Let us consider the function 


V{x,y,z,px,Py,Pz) •= Y “ A2(cos(z - a) -h 1), 


(67) 


which is a first integral of the Hamiltonian system defined by Hq, where a = aTctan{px/py). This function is used 
to estimate the distance between the invariant manifolds associated to the NHIM (see Lemma l4~2] below). Indeed, at 
every point = u^{t, x, y) G I47®(Ao) = I47“(Ao), given by (|9l), we have 

2A 


2 iA / \ 

T’{u^) = -—r “ (cos(4 arctan + tt) -|- 1) 

cosh (Ar) 4 / 


= A^ 


8 e 


2At 


-h cos (4 arctan e^”^) — 1 ) = 0 . 


( l + e 2^^)2 


Then, for every point G 147 ®(Aq) = 447“(Aq) we consider the straight line S, transversal to 417 ®(Aq) = 447“(Aq), 
given by 

S = S(u 0) = K + pV(,,p,)iP(uO) : // G R} , 

where we are using the notation V(^z,p,,)'P ■= (OjOj §^)- We denote by u® = S(n'^) n 147®(Ae) and ri“ = 

S(ri°) n VF“(Ae) the intersections of the line S with the stable and unstable manifolds of Ag, respectively. Then, there 
exist constants ^® G R and /i“ G R such that these intersections are given by 


u = { X 


= I X 


: -h Fi(r), ?/ -h Fa(r), (r) -h ,Px,Py,P^z{T) + 

' + Fi{T),y + F2(r), z°(r) -h ,px,Py,pl{T) + p'^dp^V^ 
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where := ^{u^{T,x,y)) and dp^V° = §^iu^{T,x,y)). Then, the following result states the relationship 

between V{u) and the intersections of and 14^“(Ae) for e / 0: 

Lemma 4.2. For each fixed u^, the homoclinic intersections of the stable and unstable manifolds are characterized 
by 

u^ = v}^ ^ = ^ V{u^) = V{v^). 

Proof It is clear that these implications hold from the left to the right. The converse follows from the fact that the 
function 


f{y) :=iP(uO + /iV(,,p^)iP(nO)) 


(p°(r) + 

2 


cos( 2 °(r) + ydzV^ — a) — , 


has no critical points if y is small enough. Indeed, an easy computation shows that the derivative 


~ -27—7 + sin^ (4 arctan + vr) 

cosh (At) 

does not vanish in the region {px > 0,py > 0} because A = □ 

Before proving Proposition 14. 1 1 we summarize some basic asymptotic properties of the flows and We recall 
that the dynamics of the unperturbed problem has a phase-shift u± = (a:±) y±) 0,0, 0,0), where x±, and y± are 
given by Eq. (fT^ . The trajectories on the invariant torus Tp^^py and the trajectories on the whiskers W^{Tp^^py) = 
converge exponentially to each other, with rate A as f —>• ±oo. More precisely 


- (|P^{u* +u+)\<Cie t>0, 

10°(n°) — 0°(u* + u_)| < , f < 0, 


for some constant Ci > 0. These estimates are obtained using the explicit expressions computed in Section 12.11 and 
they just reflect the normal hyperbolicity of the NHIM. Analogous expressions hold for the perturbed system. In this 
case, given u® G VE®(Ae) and u'^ G iy“(Ae), there exist points on the NHIM, u®*, G A^, that are e-close 

to their unperturbed counterparts. These points satisfy 

I0f(u") - + U^)| < , 

|0f (u“) - filin'^* + u" )| < , 


f > 0, 
f < 0, 


(69) 


for some constant C 2 > 0, where A^ = A -|- 0(e). We also need to recall some estimates that allow us to compare the 
perturbed and the unperturbed flows. The following estimates, which hold for all f G M, are standard and immediate 
to obtain (for certain positive constants Cs, 6 * 4 , 6*5 and K): 

\fUul - < C3\u^ - , 

10f (u") - 0°(u°)| < C3|u“ - , 

+ O - ftiu* + ^^+)l < < Csee-^"!*! , 

|0f (u"* + u") - 0°(u* + n_)| < C 4 |u'^* - . 

These expressions state that for e > 0 there may be unstable motions inside the NHIM and we cannot have a global 
control on the dynamics for all time. Nevertheless, we have the bounds 

< CjePf for \t\ < Celog{l/e), (71) 


with Cq > 0 , Cz > 0 and 0 < pi < 1, which is enough for our purposes. 
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Proof of Proposition \4.1\ To monitor the evolution of the function V, given by Eq. (l67l) . along the perturbed flow, we 
use the formula 

^{V{(j)l{u)) = = e{V,h}{(j)l{u)), 

where we have used that {V, Hq} = 0. Integrating this equation we obtain 

/ti 

Using (1721) with ^2 = 0, fi = oo and u = rr®, we have 

rOQ 

V{un = V{ct>Uoo{un)-e / {V,h}my^))da. 

Jo 

and using (l72l) with f 2 = 0, U = oo and u = rr®* + we get 

poo 

V{u^* + u\) = Vic^Uooiu^* + O) - W + ul))da . 

Jo 

Subtracting these expressions we obtain 

r{u^) _ r{u^* + = vicpu^iu^)) - vicpu^iu^* + <)) 

-e r [{r,h}mul) - {V,h}mu^*+u\)))da. 


(72) 


(73) 


(74) 


Now, we observe that V{u^* + rx+) = V{u* + ri+ + 0{s)) = O(e^), since both V and W vanish on the unperturbed 
NHIM Aq. Moreover, from the asymptotic properties (l69l ) it follows that 

- V{<pu^iu^* + 0)1 < C8lOoc(^^) - + OO 0. 

To study the integral term, we recall that we cannot control the dynamics on the NHIM for all time, so we consider 


r ({V,hm{nl) - {V,hm{n^*+u\)))da 

7C6log(l/e) ^ ^ 


< f Cg\(t>^^{u^) - 0^(rx®* + rx+)|dcT 

7c6log(l/e) 


< CgC 2 r e 

Jcelogil/e) 


~^'^^da = ^9^*2 ^-AeCelog(i/£) _ n(p-P^^ 


A. 


= 0{eP 


for certain constant p 2 > 0. Notice that we have used Eq. (1691 ) to derive the second inequality. Accordingly, these 
estimates and Eq. (1741) imply that 

/■Celogil/e) 

P(rx®) = -e I ({P, h}{4>Uul) - {V, + u%)))da + 0{e^) + . 

Now, we can control the quantities (j)%{u^) — and + ri+) — + ^+) using Eqs. (ITOl) and (ItTI) . so 

we write 

fC6log(l/e) 


V{un= [{V,h}{cl,^,{u^))-{V,h}{^l{u*+u+))ya + I + 0{E^)+0{e^+P^), 


where 


I <eCg 


^Celogll/e) 

(l0^(u^)-</>°(u°)| + |<^^(u®*+O-0°(rx*+rx+)|)da < 2CgCjCge^+P^ log(l/£) = 0 ( 0 ^^) 



A. Luque and D. Peralta-Salas 


27 


for certain constant 0 < /03 < 1. We conclude that 

V{u^) = -ej^ - {P,h]{<tP,{u*+u+)))da + 0{e^+P) , 

for some constant p > 0. Here we have used the bound 

r {{V,h}{fyu'‘)) - {V, + «+)))d<T = 0(e»). 

7c6log(l/e) ^ ^ 

Finally, obtaining a similar formula for P(n“) and subtracting, we obtain 

V{u^)-V{u^) = e P (^{V,h}{cP^,{u^))-{V,h}{^l{u*+u±))ya + 0{s^+n. (75) 

Recalling that the unperturbed flow satisfies Eq. (fTTl) . we can write 

= dxhP Fi{t + (j) + dyh^ F 2 {t + ct) + d^hP iP{T + a) + dp^hPP z{t + cr), 

where we are using the notation := ||((/)° (rt*^)). Now we observe that -Fi(r) = — sin(a) — sin( 2 ;°(r)), and 
F 2 {t) = — cos{a) — cos(z°(r)), so we obtain 

^(/i((/>° (tt°))) = - [sin(a) + sin( 2 ;°(T + a))]dxh^ - [cos(a) + cos(z°(r + a))]dyh° 

+ P°(t + a)dzhP + (pa;Cos(z°(r + a)) -py sin(z°(r + cr)))dp^hP . 

Using the definition of V in (1^ . we end up with 

Hence, the expression (1751) is equivalent to 

P(u") -V{u^) = -e-^L{T,x,y,px,Py) +0(e^+^), 

where we have considered the expansion h{u,e) = h{u,0) + 0{s) and used the definition of L in Eq. (1^ . By 
Eemma |4^ homoclinic intersections are characterized by the condition V{u^) = V{vF). Therefore, we conclude 
that the existence of homoclinic intersections is given, at first order perturbation theory, by the zeros of a directional 
derivative of the Poincare function L(t, x, y,Px,Py), as we wanted to prove. □ 


Eet us consider the Hamiltonian of the ABC system, written as = Hq + eHi + e^H 2 , where Hq, Hi and H 2 
are given by Eqs. (fTTl) - (fT^ . Then, using expressions (fTOl) and (fTTl) . the Poincare function L has the form 

L{t, X, y,Px,Py) = Ml cos X + M 2 cos y + M 3 sinx + M 4 siny , 

where the coefficients Mj = Mi{T,Px,Py) are given by the integrals 


Ml := B 



M2 := 




cos z*) sin(x± + ujia) — {py — cos z^) sin(Fi + uiia) — p^ cos(Fi + uiia)jda , 
sinz*) cos(y± + uj2cr) - {px - sin2°) cos(F2 + cc2(t) - p° sin(F2 + W2 Ct) idcr, 


(76) 

(77) 


M3 := B 
Ml := C 




cos z*) cos(x± + LOicr) — (py — cos z^) cos(Fi + tciu) + p^ sin(Fi + wicr) J da , 
sin z^) sin(F2 + uj2a) — {px — sin 2*) sin(y± + uj2a) — p^ cos(T2 + 0J2O') ) da . 


(78) 

(79) 


Here Fi = Fi{t + a) and F 2 = F 2 (t + a) are given by Eq. (fT^ . and z^ = z^(t + a) and p° = Pzi'^ + are given 
by Eq. dh]). 
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4.2 The Scattering map 

It is convenient to introduce the notation 


C{x - U}lT,y - UJ2r,Px,Py) ■■= UJlT,y - UJ2T,Px,Py) ■ 


(80) 


Since the properties of the unperturbed flow imply that 

L{0,X - UJiT,y - U}2T,Px,Py) = L{t, X, y, Px, Py) , 
we can consider the critical points of the function 

Tl- > C{X - UJlT,y - UJ2T,Px,Py) (81) 

in order to study the homoclinic intersections. 

Then, we introduce the domain 2? C x X C x in hypothesis A 2 of Theorem 12.21 such that for each 
{x,y,Px,Py) inV, there exists a unique critical point t* = t* {x,y,px,Py) of the map (1^ defining a smooth function 
on T>. This implies that the points 

{x,y,z^{T*),px,Py,pl{T*)) + 0{e) e iy"(A£) rti ^"(Ae) 

define a manifold Tg, called homoclinic manifold. The scattering map associated to Tg is defined in a domain C T> 
in the following way (see ifT^fTTl l: 


Se : p£,b cv ^ C T2 X 

Ub I— Uf, 

with Uf = Se(rtb) if and only if there exists u € Tg such that 

\(Ptiu) - -^0, f-^OO, 

!</’!(«)- I -^ 0 , t^-OO, 

where is the parameterization of the perturbed NHIM Tg introduced in Section 13.11 Since the parameterizing 
variables (x, y,Px,Py) and the phase space variables coincide up to order e^, we can safely assume that they are the 
same. The sets b and P^ f are defined as: 

^e,b := 1 ^ {'«b}, '■= 1 ^ 

ueFe neFe 

The scattering map relates the past asymptotic trajectory of any orbit in the homoclinic manifold to its future asymp¬ 
totic behavior. 

The scattering map (1^ is exact symplectic (see ifTTl l and it is given by the time-1 flow of the Hamiltonian 
function 

5, = 5o + e5i + 0{e^) 

where S'o corresponds to the unperturbed outer dynamics, and Si is given by the Poincare function (1^ evaluated at 
T = T*. Notice that the unperturbed scattering map for the ABC system satisfies 

ui = ub + u+-u-, u± = {x±,y±,0,0,0,0), 

where x±, and y± are given by Eq. (fT^ . Hence, we obtain the following expressions for Sq and Sp. 

So = -8 {pI+pI)^^'^, Si{x,y,px,Py) = C{x - ujiT*,y - UJ2T* ,px,Py), 
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where r* = T*{x,y,Px,Py) is the critical point of the function (IMT) that has the following expression for the ABC 
system: 

r 1 —^ Ml cos(x — ioir) + cos{y — uj 2 t) + M 3 sin(x — ujit) + M 4 sin(y — uj 2 t) , 

where M? := Mi{0,px,Py) are obtained evaluating the integrals (1761) . dTT]) . (17^ . and (17^ . Finally, we discuss some 
conditions that allows us to justify that there exists a domain V where the above construction is well posed for the 
ABC system. We fix a value of {px,Py) and notice that the function 

(x, y) 1 -^ £{x, y) = M® cos(x) + cos(y) + M 3 sin(x) + M 4 sin(y), 

has four critical points (xc, yc) given by 


M3O 

Xc = arctan 


yc = arctan —^ . 

^ MO 


(83) 


It is easy to check that these critical points are nondegenerate provided that M 4 and M 3 do not vanish simultaneously, 
and the same for M^ and M 4 . This implies, in particular, that 13 ^ 9 and (7 / 0, as required in the statement of 
Theorem 12.21 Hence, we observe that we are in the same situation considered in EOll . where the existence of r* 
was justified in detail using the tangential intersection of straight lines in the direction { 001 , 002 ) with the regular level 
curves of £(x, y), which are periodic curves which fill out a region bounded by the level curves containing the saddle 
points. 


5 Combination of inner and outer dynamics 

In this Section we conclude the proof of Theorem 12.21 and we give explicit formulas for the condition (fTSl) . To this 
end, we combine the inner and outer dynamics. In Proposition l3.9l we showed that the invariant tori (both primary and 
secondary) of the ABC system are given by the level sets of a couple of functions. This couple defines an -valued 
map thaf will be denoted as all along this section. The scattering map described in Section l4~2l transports the level 
sets of Fc onto the level sets of F^ o Then, following ifTSlfT^ . it turns out that (c.f. Lemma 10.4 in ifT^ l given two 
manifolds Si, S 2 C Ag that are invariant under the inner dynamics, if Si intersects transversally S£(S 2 ) in Ag, then 
iti VFsj- This is the main ingredient to create heteroclinic intersections between the KAM tori in A^. 

To characterize the action of the scattering map on the level sets of a given function, we follow the computations 
in in (which are also used in ITSlITOl l. Given a function F = Fq + eFi + £^F 2 + ... we can approximate F o as 

Fosc = F + {F, Sc} + 0{e^) = F + {Fq + eFi, 5o + eSi} + 0{e^) 

= F + {Fo, So} + e({Fi, Sol + {Fo, Si}) + 0{e^). (84) 

It is worth mentioning that this expression does not correspond with the expression obtained in mmiTOl . due to 
the presence of a phase-shift in the unperturbed problem. We observe that Fq and So depend only on the momenta, so 
we have {Fo, So} = 0. 

Let us consider a function F : x —>• that defines the level sets 


— Fifi -t- eFi^i + 0{e‘^) — ei, (85) 

F 2 = F 2 P 4- eF2^i -|- 0 {e‘^) = 62 , (86) 

and the transformation of these level sets by means of the scattering map 

Fiosc = Fi,o + eFi,i + e({Fi,i, So} + {Fi,o, Si}) + 0{e^) = e[ , (87) 

F2 0 Sc = F2,o + sF2,i + e({F 2 ,i, So} + {F 2 ,o, Si}) + 0{e^) = e '2 , ( 88 ) 
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for certain e'l, 62 - We will use Eqs. (IMl) - ([88]) to determine transversal intersections between these level sets. Indeed, 
if we subtract these expressions, we have 


{Fi,i,5o} + {Fi,o,5i} + 0(e) = ^^^, (89) 

{F2,i,So} + {F2,o,Si} + Oi£) = ^^^. (90) 

Then, if we use Eqs. (l85l) and (IMl) to write px = Px{x, y, ei, 62 ) and py = Py{x, y, ei, 62 ), and we introduce these 

expressions into Eqs. (I89l ) and (l90l ). it turns out that we will have intersection as long as ^—— and ^—— are small 

enough, close to the non-degenerate solutions of 

{Fi,i,5o} + {Fi,o,5i} = 0, (91) 

{F2,i,So} + {F2,o,Si} = 0. (92) 


The non-degeneracy condition, which implies that the intersection is transversal, reads as 



{Fi,i 

5o} + {Fi,o,5i}) 

U 

{Fi,i,5o} + {Fi,o,5i})\ 

v^( 

{^2,1 

So} + {F2,o,Si}) 

u 

{F2,1,So} + {F2P,S,}) j 


(93) 


for each point at the intersection of the level sets of F. Eet us remark that the condition (l9^ is evaluated by fixing px 
and Py by means of Fi^o = ei and F 2 fi = 62 , where we have neglected the 0(e) terms. 

Remark 5.1, It is easy to check that, in the non-resonant case, the matrix in Eq. (193! ) is symmetric. This is a con¬ 
sequence of the geometric structure of the problem, since the functions Fi^i and ^ 2,1 obtained by means of the 
change of variables (IMl)- 


Einally, let us express the condition (l9^ in a explicit way for the case of the ABC system. We consider separately 
the non-resonant and the resonant zones: 


• Non-resonant region: Erom Proposition 13. 9l it follows that we have to consider the level sets (l62l) . In order to 
check the condition (1931) we introduce Fi = Fi^ + eFip + 0(e‘^) and F 2 = F 2 fi -|- £^ 2,1 + where 

Fyo ■■=Px, 

Fi 1 := — B— sinx , 

F2,0 ■=Py, 

F2,i := -C^cosy. 

A direct computation shows that 

{Fi,i,S'o} = 4 :B^Px{pI+pI)~^^'^cosx , 

{Fi,o, 5 'i} = - (1 - cos(x - uJiT*) MIuj 2 ^ cos{y - uj 2 T*) 

+ (1 - sin(x - uiiT*) - M^uj2^ sin(y - uj2T*) , 

{F2,i,So} = - 4 C^py(pl + pI)-^/^ siny , 

{F2,0, 5'i} = C0S(X - tClT*) - m 2 (l- W2^) COs(y - UJ2T*) 

- sin(x - uJiT*) -h m2 (^1 - ^ 2 ^) sin(y - uj2T*) , 
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and a straightforward but cumbersome computation allows us to compute the 2x2 matrix in Eq. (l9^ . which 
reads as 

/ Ai A2 \ 

V As As J ’ 

where the coefficients have the expressions 

Ai:=i:{{F,,i,So} + {Fs,o,Si}) (94) 

= - +p 2 )- 3/4 [MgVr*^ + Mi°(l - uiiT*f] cos{x - UJiT*) 

+ [M 2 uJ 2 T*^ + M^{u! 2 T*f] COs{y - UJ2T*) + [-M%iT*^ + M ^{1 - OJlT*f] sin(x - UJlT*) 

+ [-M^uj2tI^ + Ml{uj2Tlf] sm(y - wsr*), 

As := ^({Fi,i,5o} + {Fi,o,5i}) = ^({Fs,i,5o} + {Fs,o,Si}) (95) 

= [Miu:iT*y - M^Vr*(1 - ujit*)] cos(x - ujit*) 

+ [Mluj2r*y - M^uj 2T*{1 - uj2ry)] cos{y - wst*) 

+ [—AfiWiT*j^ — MsUJiTy{l — cciT*)] sin(x — tcir*) 

+ [-M2U}2T*y - M^uj 2T*{1 - uj2Ty)] sm{y - Uj2r*), 

A3:=^({F2,i,So} + {Fs,o,Si}) (96) 

= “ ^C^pyipl + cosy + [MsUJiTyy + M°(LUir*)^] cos(a: - ujit*) 

+ [M^uj 2 Tyy + M^(l - wsr*)^] cos(y - wsr*) + [-M'^uJiTyy + M^(a;iT*)^] sin(x - ccir*) 

+ [-M2U}2Tyy + M^(l - Uj2ryf] sm(y - wst*) . 

Here the subscripts in r* denote, as usual, partial differentiation. Then, the transversality condition in the 
non-resonant region, using the functions A,, reads as 

AiAs - A2 / 0. (97) 

• Resonant region: From Proposition l3.9l it follows that we have to consider the level sets (l6^ . In order to check 
the condition in Eq. (|9^ we introduce Fi = Fi^o + and Fs = Fs^o + £-^ 2,1 + Cl(e^), where 

i^l ,0 ■=Px+Py, 

Fi,i := - sin X - cosy, 

F2,o :=| + :^ + (l + V2er')(f-p,)', 

Fs,i := 217^(1 + V2er') (f - Py) cos y . 

Then, the matrix in Eq. (|93] ) has the form 

f Ai As \ 

V As A4 y ’ 

with 

Ai := i({Fi,i,So} + {Fi,o,Si}) = -4F^p,(p2 +p2)-3/4g.^^ 

+ [Mluiir*^ + T*y) + M {’(1 - a;iT,^)(l - wi(r* + t*))] cos{x - ujit*) 

+ [M4°a;s(T,^,^ + T*y) - M^uj2T*il - uj2{t* + T*))] cos(y - WsT*) 

+ [-A/iWi(r*,, + T^y) + M^{1 - uJiT*){l - cji(r* + t*))] sin(x - ujit*) 

+ [-A/s°ws(r*,^ + T*y) - M^uj2T*{1 - a;s(r* + t*))] sin(y - tcsr*), 


(98) 
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A 2 := + {Fi,o,Si}) = -4C^^py{pl+pl)-^/\osy (99) 

+ [MluJl{T*y + Tyy) - M^UJlT^il - W 1 (T^ + Ty ))] COs{x - WlT^) 

+ [M2u}2iT*y + T*y) + M^{1 - uj 2 T *){1 - 0 J 2 {t* + T*))] cos{y - U}2T*) 

+ [-M^UJi{T*y + Tyy) - UJiTy - W 1 (T^ + T))] sin(x - WlT^) 

+ [-M^UJ2{T*y + Tyy) + M^(l - UJ2Ty){l - UJ 2 {t^ + Ty))] sm{y - UJ2T*) , 

A 3 := £i{F 2 ,i,So} + {F 2 ,o, Si}) = 7 [-MgVr^ + MiVt;(1 - uit^)] cos(x - loit*) (100) 
+ 7 [-M 4 Ui 2 r*y + M^uj 2 T *{1 - W 2 T*)] cos(y - uj 2 T*) 

+ 7 [MiU!lT*y + M^COlTy{l - UIlT*)] sm(3: - UIlT*) 

+ 7 [M 2 ^ 2 T*y + M^UJ 2 T *{1 - Ul 2 Ty)] Sm{y - U 2 T*) , 

A 4 := ^{{F2,1,So} + {F2,0,Si}) = SCpyipl +p2)-3/4^(i + ( 101 ) 

+ 7 [-M^UJlTyy - M^iUlTy)'^] COS{x ■ 1 T ^ ) 

+ 7 [-M 4 u; 2 Tyy - M^(l - W 2 r*) 2 ] cos(y - a; 2 r*) 

+ 7 - M;^(a;ir*) 2 ] sin(x - luit*) 

+ 7 [M2(^2Tyy - M^il - W2T*)^] sm{y - U2T*) . 

Here we are using the notation 7 ;= 2(1 + \/2) . Then, the transversal!ty condition in the resonant 

region takes the form 

A 1 A 4 -A 2 A 3 /O. ( 102 ) 

Putting together the information gathered on the inner and the outer dynamics, we can construct chains of invariant 
tori giving rise to large motions in the action space. Assume that V is the domain introduced in Section 14.21 It is 
obvious that we can safely assume, by shrinking V if necessary, that the domain b of the scattering map (1^ 
coincides with V. Then, we assume that we can choose a constant L such that such that for every (x, y,Px,Py) G F 
we have 

r A 1 A 3 -AI/O, if \Px-Py\>L, 

\ A 1 A 4 - A 2 A 3 / 0 , if \px -Py\ <L. 

This condition is precisely Eq. ([TSll in Theorem l2.2l It then follows that we can find a sequence which 

are at a distance 0 {e) from each other and that satisfy Se{Ti) (ti Ti+i- By applying Lemma 10.4 in ifT^ . it turns out 
that these tori satisfy Wp. iti , that is, they form a transition chain. The claim of Theorem |2!2] (Arnold diffusion) 
then follows from the general theory presented in ifT^fT^ . 

6 Rigorous verification of the hypotheses of the main theorem 

In this section we illustrate the effective verification of the hypotheses of Theorem l2.21 thus obtaining Arnold diffusion 
in the ABC Hamiltonian system (fT4l) . We rigorously evaluate the involved functions and obtain rigorous bounds for 
the critical points with the help of the computer. Our approach is as simple as possible, in the sense that we do not 
pretend to present a fast and efficient methodology to study large regions of the phase-space systematically. Our 
interest here is to convince the general reader that the hypotheses of Theorem I2.2l can be rigorously checked with the 
help of a computer. 

In rigorous computations, real numbers are substituted by intervals whose extrema are computer representable real 
numbers. That is, when implementing interval operations in a computer, the result of an operation with intervals is an 
interval that includes the result. The reader can consult the recent introductory book ll44l on rigorous computations. 
All the computations presented in this section have been performed using FILIB |[^ that uses double precision 
arithmetics. 
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Rigorous bounds of the Melnikov coefficients are obtained in Section [ 6 A] In Section we control the critical 
point T* and its derivatives with respect to the angles {x, y). Finally, in Section [631 we present a direct application of 
the previous ideas and we describe the implementation details giving rise to Corollary 11.11 


6.1 On the evaluation of the Melnikov coefficients 


Given certain values of px and py (represented using interval arithmetics), we are interested in the rigorous evaluation 
of the Melnikov coefficients Mi, M 2 , M 3 and M 4 in Eqs. (17^ - (|79]). Let us recall that we are particularly interested 
in the values Mf = Mi{0,px,Py) in order to check the hypotheses of Theorem 12.21 If we denote by fi{a) the 
function that we have to integrate to evaluate M?, and we introduce the notation /j _|_(cj) = /i(cr) for cr > 0 and 
fi-{a) = fi{a) for a < 0 , we can write the expressions for M° as follows 

/ oo p—a /*0 pa poo 

fi{a)da = / fi-{cr)da + / fi-{a)da+ / fi^+{a)da + / fi^+{a)da, 

■OO J—00 J—a Jo Ja 

where a > 0 is a constant that will be fixed lafer. The infegrals af infinify (called tails from now on) will be bounded 
using fhe asympfofic properties discussed in Section |2T] Of course, one can obfain general formulas for fhe fails in 
terms of a uniform confrol on fhe Hamiltonian and fhe Lyapunov exponenf. However, in fhis case, we present specific 
formulas for fhe ABC system giving rise fo sharper esfimafes of fhe fails. This allows us fo keep fhe modulus of fhe 
fails under a prefixed folerance using a small value of a. 


Lemma 6.1. The following bounds hold for the ABC system: 

rO pa 

fO 


Mi- fi-{a)da- /i,+ (cr)dcr 
J—a J 0 


< 5],, 


where Sj are given by 


51 = S 3 := B^\py - cosz* 

4 

5 2 = S 4 := C^\px - sinz* 


log(l + e ^^“)sinz 


-2( 


arctan 


-.—Xa 


- 1 


g-Aa _j_ I 
—Xa _ ^ 


TT 

+ — ) cos z’ 


( / 0 Xa _ \ ^ \ 

arctan _xa -p 1 j 


+ -BSq , 

+ CTjQ , 


+ 4i?| arctan(sinh(Aa)) — vrl. 


whereto — 8^ + 3^ 

Proof We will only consider fhe case of because fhe confrol of fhe fails of M^, Mg and M 4 is complefely 
analogous. We firsf splif fhe function /i _|_ info fhree ferms as 

f _ f(t) I f(2) I f{3) 

In fhis splitting each term is given by 

:= B{py — cos2*)(sin(x+ + uiia) — sin(Fi((T) + uiia)) , 

:= .Bsin(Fi(cr) + cjicr)(cos2:'^(cr) — cosz*), 
ffX ■■= - Bp°{cr) cos(Fi((t) + uia). 

Then, a slraighlforward compulation shows lhal 


poo 

J a 

poo 

J a 

poo 

/ fiJ{(y)d(j 

J a 


< B\Py — COS z* 


< B 


< B 


poo 

/ (4^ 

J a 

f 


POO 

/ {Fi{a)-x+)da 
J a 


arctan — 2Tr)da 


2A 


cosh(Acr) 


da 


(103) 

(104) 

(105) 
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By the asymptotic properties discussed in Section 12.11 we know that these three integrals are convergent. Next, we 
give some explicit expressions to control the above integrals. To this end, we use the expression of the primitives of 
the functions that we are integrating. First, we introduce 


g{a) := j ( 2 ^ 1 ( 17 ) — x+)da = ^ log(cosh(Acr)) sin z* — arctan(tanh(^))-cos 2 * — x+cr , 

which allows us to control Eq. (11031 ) in terms of the expression g{oo) — g{a). However, the direct evaluation of this 
expression with a computer presents a huge rounding error. A more suitable formula is obtained using the limit 


lim logfcoshlAa)) = lim log(- 

(J—>CX3 C7 —>00 

which allows us to control the term (11031) as follows 


-) = lim (Afj — log 2), 

(T—)-(20 


f 


fi^l{cr)da 


< Bj^\py - COSz* 


log(l + e ^^“)sinz 


-2( 




arctan 


Aa 


- 1 


:i —Ao. 


+ 1 


TT 

+ — I COS z’ 


Using Taylor series, it is easy to check that the term (11041) is bounded as 

3 —Aa ^ g—3Aa 


f 

J a 


(4 arctan — 2'K)da 


< 4- 


A 


+ 


3 A 


and that the term (11051) is estimated as 



2A 

cosh(AiT) 


da 


< 2| arctan(sinh(Aa)) — tt 


Analogously, we can estimate the term 


thus proving the lemma. 


f fi_{a)da ^ 
J — 00 


□ 


We use Lemma lhdl to rigorously control the Melnikov coefficients. Specifically, we evaluate directly the obtained 
expressions of Sj using interval arithmetics. The integrals and /*_+ are controlled using Simpson’s rule 

with rigorous bounds on the error, which are obtained using explicit formulas for the 4th-order derivatives of the 
functions fi(a) computed with a symbolic manipulator. 

Next, we illustrate the rigorous evaluation of the Hypothesis Ai of Theorem 12.21 Obviously, it is enough to 
consider the coefficients for the parameters B = C = 1 and store conveniently the obtained values. If we are 
interested in other values of B and C we simply have to scale the previously computed values. In Table[T]we present 
some rigorous enclosures of the coefficients corresponding io B = C = 1, px & [0.4,0.4001] and different 
interval values of py. To control the tails we use Lemma lOl with a = 20 and to enclose the finite integrals we use 
Simpson’s rule with 130 subintervals. This implementation parameters are enough to guarantee that the coefficients 
Mf do not vanish for a non-empty set Z of momenta. We can obtain a similar result for a much larger domain Z by 
systematically performing this computation. 

Finally, let us illustrate how to check Hypothesis A 2 of Theorem 12.21 It turns out that, for px £ [0.4,0.4001] and 
all the intervals py in Table[T] the four critical points of the map (x, y) 1 —)• C{x, y) in Eq. ([83]) are non-degenerate and 
given by two maxima and two saddle points, so that we can use the arguments in ll20l to see that there is a unique 
smooth critical point r* of the function r 1 —)• C{x — uoiT,y — L 02 T,px,Py)- In the following section we discuss the 
rigorous enclosure of r* and its derivatives. 
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Py 

MV 

MV 

MV 

MV 

[0.1,0.1001] 

[-5.1237,-4.9193] 

[-11.314,-10.972] 

[-1.282,-1.0754] 

[1.4611,1.807] 

[0.2,0.2001] 

[-5.3464,-5.1279] 

[-9.185,-8.8668] 

[-2.2718,-2.0471] 

[-3.0694, -2.7531] 

[0.3,0.3001] 

[-5.911,-5.6582] 

[-6.5333, -6.2423] 

[-3.3557, -3.0967] 

[-5.5463, -5.2556] 

[0.4,0.4001] 

[-6.4914,-6.2559] 

[-4.4201,-4.1843] 

[-4.4214,-4.1829] 

[-6.4917,-6.2557] 

[0.5,0.5001] 

[-7.0075, -6.7721] 

[-2.9704,-2.7501] 

[-5.2892, -5.0544] 

[-6.7487, -6.5249] 

[0.6,0.6001] 

[-7.4463, -7.2257] 

[-1.988,-1.7805] 

[-5.8843, -5.6657] 

[-6.701, -6.4904] 

[0.7,0.7001] 

[-7.8251,-7.6315] 

[-1.3078,-1.129] 

[-6.2495, -6.0596] 

[-6.5182,-6.3291] 

[0.8,0.8001] 

[-8.1753,-7.9784] 

[-0.85053, -0.66663] 

[-6.4659,-6.2811] 

[-6.2954, -6.0987] 

[0.9,0.9001] 

[-8.4879, -8.2889] 

[-0.53277, -0.3423] 

[-6.575,-6.3929] 

[-6.0543, -5.8485] 


Table 1: We show rigorous enclosures of the coefficients Mf, for the values B = C = 1 andpa; € [0.4, 0.4001]. 


6.2 On the evaluation of the critical points 

We discuss a simple methodology to rigorously enclose the critical points of the function 

r I—1 C{x - LOlT,y- U} 2 T,Px,Py) 

= Ml cos{x — (jOit) + M 2 cos{y — u} 2 t) + M 3 sin(a: — wir) + M 4 sm{y — lo 2 t) , (106) 

where {x,y,Px,Py) G x X, and the coefficients M? = Mi{0,px,Py) are given by Eqs. (176)) - (17^ . The critical 
points T* = T*{x,y,px,Py) of (1106b are characterized by the zeros of the function 

Q{t) := tdi ^ — Ml sin(x — loit) + Mg cos(x — tuir)^ + 0^2 ^ — M® sm{y — uj 2 t) + M 4 cos(y — W 2 t)^ , 

which, of course, depends on the variables {x,y,px,Py)- There are several techniques in computer-assisted proofs 
that allow us to study solutions of nonlinear equations as above, e.g. the interval Newton method Il4^ . However, in 
the following discussion we choose to use the simplest possible method with the aim of convincing a reader that is not 
familiar with computer-assisted methods. More advanced techniques will give the possibility to validate large regions 
of phase-space with a reduced computational cost. This is not the aim in this article, since the required techniques are 
not related with the ideas that we want to highlight and they would require to provide a larger amount of computational 
and implementation details. To enclose r* we proceed using a bisection-like procedure: 

• Given {x,y,px,Py), that may be numbers or interval values, we enclose M? = Mi{0,px,Py) following Sec¬ 
tion | 6 T] 

• Given an integer > 1, we consider an increasing sequence of interval values Tj G [(f — 1) /N,i/N]. Then we 
compute Q{Ti). While 0 ^ Q{Ti) we increase the index i, until we obtain an interval such that 0 G Q(ri). This 
gives a lower estimate for t* . 

• After the previous computations, we continue the process of increasing the index i and computing Q{Ti). When 
we reach an interval such that 0 ^ Q{Ti) then we have obtained an upper estimate for r*. 

As a result of the above procedure, we obtain an interval enclosure of the critical point r*. Then, the derivatives 
of T* with respect to {x, y) are computed using the following equations: 

Qa{T*) Qa/t(r*) + Q„.(T*)r^ + (Q^,(r*) + Q,,(r*)T|)r* 

Qr{T*) 

where the subscripts denote, as usual, partial differentiation, and a and /3 can be chosen to be x or y. For the case of 
the ABC system, we have 

Qx('r) = — Ml cos(x — ojit)uji — Mg sin(x — a;iT)a;i, 

Qyir) = - M 2 cos(y - 0 J 2 t)uj 2 - M 4 sm{y - uj 2 t)oj 2 , 










A. Luque and D. Peralta-Salas 


36 


Qrir) = cos{x — + Mg sin(x — cos{y — W2T)a;| + M4 sm{y — u:2t)oj2 , 

Qxxi^) = Ml sin(x — a;ir)a;i — Mg cos(x — ojit)uji , 

Qyyir) = M2 sm{y - UJ2t)uJ2 - M4 COs(y - 0 J 2 t)ijJ 2 , 

Qa;T('r) = — Ml sin(x — ujit)ujI — Mg cos(x — ujit)oj{ , 

Qyrir) = — M2 sin(y — uj2t)uj2 — M4 cos(y — oj2t)lo2 j 

< 3 Tr('r) = Ml sin(x — ujit)coI — Mg cos(a; — coit)uiI + M^ sin(y — u!2t)uj2 — M4 sin(y — u:2t)ijj\ . 


In Table |2] we illustrate these computations taking B = C = 1 and considering the values of px and py used in 
Table [U In all the cases, we fix the angles as x = y = 0, and we compute the critical point r* using N = 100. 


Py 

T* 


'^y 

[0.1,0.1001] 

[2.34,2.43] 

[0.74112,0.91778] 

[-0.41917,-0.3401] 

[0.2,0.2001] 

[2.71,2.81] 

[o.47369, 0.73294] 

[0.29594, 0.45414] 

[0.3,0.3001] 

[2.38,2.45] 

[0.436II, 0.54436] 

[0.41137,0.51753] 

[0.4,0.400lj 

[2.1,2.16] 

[0.40377, 0.50552] 

[0.40346,0.5051] 

[0.5,0.5001] 

[1.9,1.96] 

[0.37234, 0.4717] 

[0.39116,0.50717] 

[0.6,0.600lj 

[1.75,1.81] 

[0.33403, 0.42893] 

[o.38246, 0.51976] 

[0.7,0.700lj 

[l.61,1.67] 

[0.27007, 0.3617l] 

[o.38473, 0.54755] 

[0.8,0.8001] 

[1.48,1.55] 

[0.1825,0.2837] 

[0.37799, 0.60038] 

[0.9,0.9001] 

]l.36,1.44] 

[0.092554,0.18709] 

[0.37117,0.65404] 

Py 

T* 

' XX 

T* 

‘xy 

'^yy 

[0.1,0.1001] 

[0.45478,1.4543] 

[-0.25048, 0.00020579] 

[-0.40322, -0.24827] 

[0.2,0.2001] 

]o.48773,1.5132] 

[0.22382,0.81026] 

[-0.363, 0.057794] 

[0.3,0.3001] 

[0.41841,0.76705] 

[0.32713,0.5756] 

[0.18447,0.45382] 

[0.4,0.4001] 

[0.24271,0.51779] 

[0.34875, 0.58267] 

[o.40865, 0.70313] 

[0.5,0.5001] 

[0.066261,0.29543] 

[0.36568, 0.62845] 

[0.62318,1.0275] 

[0.6,0.6001] 

[-0.11457,0.079559] 

[o.38581, 0.69138] 

[0.85943,1.4404] 

[0.7,0.7001] 

[-0.32399,-0.12284] 

[o.37369, 0.72931] 

[1.0946,1.9373] 

[0.8,0.8001] 

[-0.53479, -0.27236] 

[0.29473,0.75013] 

[l.2426, 2.6407] 

[0.9,0.9001] 

[-0.68739, -0.39066] 

[0.18037,0.6502l] 

[1.3104, 3.4064] 


Table 2: We show rigorous enclosures of the critical point r* and its derivatives with respect to (x, y) for the values B = C = 1, 
Px € [0.4, 0.4001], and x = y = 0. 


6.3 On the verification of the transversality conditions 

With the rigorous estimates obtained in Sections 16.11 and 16.21 let us now explain how to use interval arithmetics 
to check the transversality condition in Hypothesis Ag of Theorem 12.21 This consists in enclosing the functions 
{Ai}j=i 2,3 and {Ai}i=i 2 , 3,4 given by Eqs. (l94l) - (IIOII) . In order to check the condition in the resonant region, we 
observe that A 3 and A 4 are proportional to px — Py so that they tend to zero when we approach the resonance. For 
this reason, we eliminate the factor px — Py in the computation of the expression A 1 A 4 — A 2 A 3 . This allows us to 
check that the condition holds in any tubular neighborhood that is close enough to the resonance. 

An important observation it that we have the freedom of choosing the angles (x, y) to evaluate the critical point r*. 
In fact, there is no optimal way to choose the angles (x, y) in order to verify the transversality conditions. The reason is 
that optimal values selected numerically may fail to fulfill such condifions when rigorous interval operations are used. 
This is because fhe enclosed value of Qt{t*) may be very close fo zero (or even confain fhis poinf), fhus producing a 
large enclosure in fhe evaluation of fhe condifions. Our experience in fhis problem is fhaf choosing random values of 
(x, y), unfil we reach a suifable pair, is fhe simplesf and fasl sfrafegy. 

For example, in Table [3] we presenf some rigorous enclosures of fhe fransversalify condifions in Hypofhesis A 3 of 
Theorem 12.21 corresponding to B = C = I, px ^ [0.4,0.40001], and differenf interval values of py. To confrol fhe 
fails we use Femma lhdl wifh a = 20 and fo enclose fhe finite integrals we use Simpson’s rule wifh 300 subintervals. 
To obfain fhe crifical poinf and ifs derivafives, we use fhe approach described in Section 16.21 wifh N = 300. This 
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implementation parameters are enough to guarantee that the functions A 1 A 3 — A| and Ai A 4 — A 2 A 3 do not vanish 
for a non-empty set I of momenta. By computing simultaneously the condition in the non-resonant region (4th 
column of Table |3]l and in the resonant region (5th column of Table O, it is clear that we can select a number L > 0 
that allows us to obtain diffusing orbits crossing the resonance. 


Py 

X 

y 

A 1 A 3 - A2 

A 1 A 4 —A 2 A 3 

Px-Pv 

[0.1,0.10001] 

5.2923 

0.93117 

[-26.899, -3.5905] 

[2.3401,101.72] 

[0.2,0.20001] 

2.7665 

0.55732 

[8.3643, 23.588] 

[-68.148,-6.1705] 

[0.3,0.30001] 

1.1969 

0.37322 

[26.796,48.326] 

[-135.08,-61.886] 

[0.4,0.4000lj 

3.7869 

4.19530 

[11.818,28.517] 

[-73.026, -3O.O82] 

[0.5,0.5000lj 

4.9160 

0.61701 

[-13.819,-6.4592] 

[18.387,34.8] 

[0.6,0.6000lj 

1.7342 

2.53840 

[-17.72,-3.2241] 

[8.7055, 41.964] 

[0.7,0.70001] 

4.7928 

5.74470 

[-16.15,-6.062l[ 

[13.441,36.899] 

[0.8,0.8000lj 

5.0542 

6.19650 

[-20.256,-6.1962] 

[12.115,45.704] 

[0.9,0.90001] 

1.5249 

4.18960 

[-6.1939,-0.93573] 

[0.38193,12.601] 


Table 3: We show rigorous enclosures of the transversality conditions of Theorem l2.2l for the values B = C = 1, Px d [0.4, 0.40001], for 
different interval values of py. The critical point r* and its derivatives are evaluated at the points {x, y) shown in the 2nd and the 3rd columns. 


Finally, we describe the implementation parameters of the CAP that lead to the result stated in Corollary 11.11 We 
take .6 = 10 and C = 0.1 and we divide the set Z = [0.1,0.9] x [0.5, 0.9] in subsets of size 10“'^ x 10“^. For every 
subset, we use Lemma lhdl with a = 20 and to enclose the finite integrals we use Simpson’s rule with 130 subintervals. 
To obtain the critical point and its derivatives, we use the approach described in Section [ 6 ^ with N = 100. For all 
these sets we obtain that the function Ai A 3 — A| does not vanish in Z, and the function Ai A 4 — A 2 A 3 only vanishes 
on the resonant line px = Py 
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